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Abstract. I review the current status of quasar black hole (BH) mass estimations. 
Spectroscopic methods have been developed to estimate BH mass in broad line quas- 
ars to an accuracy of ~ 0.5 dex. Despite their popularity, significant issues and con- 
fusion remain regarding these mass estimators. I provide an in-depth discussion on 
the merits and caveats of the single-epoch (SE) virial BH mass estimators, and a de- 
tailed derivation of the statistical biases of these SE mass estimates resulting from 
their errors. I show that error-induced sample biases on the order of a factor of several 
are likely present in the SE mass estimates for flux-limited, statistical quasar samples, 
and the distribution of SE masses in finite luminosity bins can be narrower than the 
nominal uncertainty of these mass estimates. I then discuss the latest applications of 
SE virial masses in quasar studies, including the early growth of supermassive black 
holes, quasar demography in the mass-luminosity plane, and the evolution of the BH- 
host scaling relations, with specific emphases on selection eff'ects and sample biases in 
the SE masses. I conclude that there is a pressing need to understand and deal with the 
errors in these BH mass estimates, and to improve these BH weighing methods with 
substantially more and better reverberation mapping data. 
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1. Introduction 



Shortly after the discovery of quasars at great cosmological distances (ISchmidtlll963h . it was 



realized that the energy required to power these luminou s and compact sources must be of grav 



itational origin rather than from nuclear reaction (e.g., iHovle & Fowler! 1 1963t ISalpeteiill964 



*In this review I use the terms "quasar" and "Active Galactic Nucleus (AGN)" interchangeably to refer to all active 
supermassive black holes, although traditionally quasai's are loosed defined as the luminous (L(,oi > lO'- Lq) subset of 
AGNs. By default I use quasai's and AGNs to refer to unobscured (Type 1), broad-Hne objects unless otherwise specified. 
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Zel'dovich & Novikov 1964 : Lvnden-Bell 19691) . The standard picture now is that mass is ac- 
creted onto a supermassive black hole (SMBH) at the center of the galaxy, and the gravitational 
energy is released during this accretion process to power quasar activity. If the SMBH grows 
mostly via this accretion process, its mass growth rate is simply: Mbh - '^L^ddi^ - £)/{ec^), 
where Ledd = IMbh = 1.26 x IO^^Mbh/M^) erg s ' is the Eddington luminosity of the BH, 
A = Lboi/^Edd is the Eddington ratio, and e is the radiative efficiency, i.e., the fraction of accreted 
rest mass energy converted into radiation. If both A and e are non-evolving, the BH mass increases 
by one e-fold on a characteristic timescale = ec^/[(l - e)Al] ^ 4.5 x 10^ yr, also known as 
the Salpeter time or e-folding time. If quasars do not radiate beyond the Eddington limit A - 1, 
the o bserved luminosity provides a lower-limit on their BH mass (e.g., IZeFdovich & Novikov^ 
1964 ) . The discovery of lu minous quasars (with Lboi k. 10^^ erg s ') at z > 6 (e.g.. Fan et alj 



2001 



Mortlock et al.ll201 ll) then suggests that SMBHs with Mbh > 10 Mq are already formed 



in the first billion year after the Big Bang. 



In the past two decades or so, there has been tremendous progress in the demographic studies 
of SMBHs in the nucl ei of near by galaxies (for recen t reviews, see, e.g. jKormendv & Richstone 
1995l:lFerrarese & Fo rd 2005; Kormendv & Holl2013h . It has come to the consensus that SMBHs 
with masses of ~ 10^ - 10'*^ M© are almost ubiquitous at the center of massive galaxies with 
a significant spheroidal (bulge) component, and also exist in at least some low-mass galaxies. 
More remarkably, the mass of the nuclear BH is tightly correlated with the properties of the bulge 
in the local samples (e.g., Gebhardt et al. 2000; Ferrarese & Merritt 2000; Graham et al. 2001] 
Tremaine et al.ll2002tlM^coni & Hunt2003: Aller & Richstone 2007; Gultekin et al.. 2009). al- 
lowing an estimate of the local SMBH mass function by convolutions with galaxy bulge distribu- 
tion functions. These BH-bulge scaling relations promoted the notion of BH-galaxy co-evolution, 
during which the energy release from the accreting SMBH self -regulate s its growth, and impacts 
the formation and evolution of the bulge via feedback processes (e.g. JSilk & Reeslll998t iKin^ 
20031: bi Matteo et allEoOSl) . Such feedback from active SMBHs (i.e., AGN feedback) has also 
been invoked in most present-day theoretical modeling of galaxy formation, to bring better agree- 
ment with the observed statistics of massive galaxies. However, the significance of AGN feedback 
and BH-host co-evolution is still under some debate and is an active area of research. 



An elegant argu ment tying the local relic SMBH population to the past active population is 
the Soltan argument ('Soltan' 1982): if SMBHs grow mainly through a luminous (or obscured) 
quasar phase, then the accreted luminosity density of quasars to z = 0, p.,acc, should equal the 
local relic BH mass density p.: 



Jo dz Jo 



(1 - e)L 



<1)(L, z)dL K p. 



(1) 



where <^{L,z) is the bolometric luminosity function (LF) per L interval. Given the observed 
quasar luminosity function, a reasonably good match between p..acc and p, can be achieved if the 
average radiative efficiency e ~ 0.1 (e.g., Yu & Tremaine,2002 : Shankar et al.l2004 ; Marconi et al. 
20041 also see Salucci et al. 1999; Fabian 1999; Elvis et al. 2002), consistent with the m ean e 
value constrained from individual quasars with spectral fitting (e.g.. IPavis & Laoi 2011). The 
Soltan argument and its variants have been used extensively in recent years to model the growth of 
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SMBHs with constraints from the demo graphies of local BH relics and the past AGN population 
(for a recent review, see IShankarlE009l . These exercises are mainly facilitated by the advent of 
modern large-scale, multiwavelength sky surveys, which have provided large and homogeneous 
data sets many folds more than what was available twenty years ago, as well as measurements of 
the abundance and clustering properties of quasars with unprecedented precision. 

The growth of SMBHs is among the key science topics in modern gal axy formation studies 
(for a relatively complete summary of recent progress on this topic, see, e.g. . [Alexander & Hickox 
2OI2I and references therein). As one of the few fundamental quantities describing a BH, the mass 
of quasars is of paramount importance to essentially all quasar-related science: the evolution and 
phenomenology of quasars, accretion physics, the relations and interplays between SMBHs and 
their host galaxies. 



In this review I discuss the current status of quasar BH mass estimations and how these 
developments can further our understandings of the physics and evolution of SMBHs. I presume 
the reader has a basic understanding of AGNs and I will s kip elaborations on the usua l AGN 
terminologies, which can be found in AGN textbooks re.g.. lPetersonlll997 : lKroliklll999h . This 
review is mostly pragmatic without going into the detailed and sometimes poorly understood 
physics behind observations; some further readings can be found in the quoted references. 



T here are sev eral recent reviews on measuring active and ina ctive BH masses (e .g.,'Petersonl 

2010l:ICzemv & Nikolaiuk.2010; Vestergaai-d et al. 201 1: Mai-ziani & Sulenticl2012l) . which sum- 
marized some general concepts and practical procedures in measuring SMBH masses. While 
some of the common materials are also covered in the current review for completeness, the scope 
and focus of this review are different: after an introduction on BH mass measurements in ^ I 
describe in detail the caveats and statistical biases of the most frequently used BH mass estimat- 
ors in ^ in light of recent work invoking statistical quasar samples; several applications of these 
BH mass estimates to quasar studies are discussed in ^ and I conclude this review in ^with a 
discussion on future perspectives of improving BH weighing methods. A flat ACDM cosmology 
is adopted throughout this review, with Oa = 0.7, Qq = 0.3 and //o = 70 km s 'Mpc~'. 



2. Methods to Measure Quasar BH Masses 



2.1 Virial BH Masses: From Reverberation Mapping to Single-Epoch Metliods 



Reverberation mapping The broad emission line r egions in AGNs are powered by photoioniz- 
ations from the central source (e.g.. |Petersonlll997l) . an assumption now widely accepted based 



on observations of correlated broad line and con tinuum variations. In fact, this photoioniza- 
tion assumption led to the first suggestions (e.g., iBahcall et al. I ll972h of lagged broad line re- 
sponses to continuum variations, where the lag reflects the light travel time from the ionizing 
source to the broad line region (BLR). In the 1970s there were already reported lag measure- 
ments between broad emission line and continuum variations in several local Seyfert nuclei (e.g., 
Lvutvi & Cherepashchukll 19741) . Despite the low data quality that may impact the reliability of 
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the detection, these studies were among the first attempts to directly measure BLR sizes. This 
idea was later developed in greater detail by Blandford & McKee ( 1982), who suggested that by 
mapping the response function of the broad emission line to continuum variations one can in 
principle reconstruct the structure and kinematics of the BLR, a technique they coined "reverber- 
ation mapping" (RM). Today RM has become a practical and powerful tool to study BLRs (see 
reviews by, e. g.. IPetersonI 1 9931 : iNetz er & Pete rsonll997[lHorne et al.l2004l) . whose spatial extent 
(~ sub-pc) is too small to be resolved by current instrumentation. There are no w several dozens of 
AGNs and quasars (most are at z < 0.3) with average lag measurements (e.g., iKaspi et al. I I2OOOI: 



Peterson et al.ll2004t iBentz et al . 2009b), although only a handfu l of them have de cent velocity- 



resolved delay maps (e.g., Dennev et al.il2009a ; Bentz et aDl2010l : Grier et al]l2013 ) to utilize the 
full power of RM. 

RM lag measurements provide an estimate of the typical size of the BLR. If we further 
assume that the BLR is virialized and the motion of the emitting clouds is domina ted by the 



gravitational field o f the central BH, then the mass of the BH is determined by (e.g., 1999 



Wandel et al.lll999h : 

Mrm = ^=/^, (2) 

where Vvii is the virial velocity and R is the BLR size. In practice we use the width of the broad 
line, W, as an indicator of the virial velocity, assuming that the broad line is Doppler broadened 
by the virial motion of the emitting gas. The product W^RjG is called the virial product. There 
are two commonly used line width definitions, the "full-width-at-half-maximum" (FWHM), and 
the line dispersion crwae (i.e., the second moment of the line, Peterson et al. 2004). The pros 
and cons of both definitions will be discussed later. In computing the RM BH masses, both line 
widths are measured from the rms spectra from the monitoring period, thus only the variable part 
of the line contributes to the line width calculation. 

To account for our ignorance of the structure and geometry of the BLR which determine 
the relation between the virial velocity and the line-of-sight (LOS) velocity inferred from W, 
we have introduced a virial coefficient (or geometrical factor), ^ in Eqn. (|2]l. This is a big 
simplification, because the BLR structure and viewing angle determine the entire line profile, 
and the line width, being only one characteristic of the line profile, can not fully describe the 
underlying kinematic structure. Similarly it is an approximation to describe the BLR with a 
single radius R. Nevertheless, given the difficulties and ambiguities of modeling the line profile 
directly, Eqn. ^ involving line widths and / is used almost universally. For BLR cl ouds in 



randomly orientated orbits, an often quoted value is / a: 3/4 (3) if W = FWHM (crune) (INetzei 
[1990), although such a / value is derived under some simplifications and approximations and 
is not a rigorous analytic result. In practice, the value of / is now empirically determined by 
requiring that the derived RM masses are consistent with those predicted from t he BH mass- 
bulge stellar velocity dispersion (Mbh - o",) relation of local inactive galaxies (e.g.. lOnken et al 



2004 : / ~ 1.4 (5.5) for W = FWHM (crii„e). This / value is of course then the averaged value 



*So me studies define the virial coefficient differently, i.e., Vyii = /W (e.g.. iMcLure & Jarvisll20d2L iDecarli et all 
120083) . 
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for the subset of RM AGNs with bulge stellar velocity dispersion measurements. The uncertainty 
in / and the simplification of it as a single constant remain one of the major uncertainties in RM 
mass determinations, as further discussed in ^3.1.21 

The uncer tainty of the RM masses is typically a factor of a few, or ~ 0.4 - 0.5 dex (e.g., 
PetersonI 2010l) . based on comparisons between RM masses and predictions from the Mbh - cr* 



relation, and accounting for other potential systematics. 

The R - L relation Perhaps the most remarkable finding of RM observations is a tight correla- 
tion between the measured BLR size and the adjacent optical continuum luminosity Lopt (usually 
measured at restframe 5 100 A), R cc L", over ~ 4 orders of magnitud e in luminosity. This is 



known as the BLR size-luminosity relation, or the R - L relation (e.g., Kaspi et al. 200d 2005 
Bentz et al l2009ah . To first order all broad line quasars have similar spectral energy distributions 
(SEDs) from X-ray to opticaj^ so Lopt is proportional to the ionizing continuum Lion. The ioniza- 
tion parameter in a photoionized medium is = Q(H)l{Anr^cne), where Q(H) oc L is the number 
of ionizing photons from the central source per second, c is the speed of light, and tie is electron 
density. Thus a slope of a = 0.5 in the R - L relation is expected, if U and the electron density 
are more or less constant in BLRs. Alternatively, a slope of a = 0.5 is also predicted if the BLR 
size is set by dust sublimation (e.g.. lNetzer & Laoilll993h . 



Early RM work reported a slope of a ~ 0.7 (e.g., iKaspi et al. I I2OOOI) . Later work which 



carefully accounted for host starlight contamination to Lopt reported a * 0.5 (e.g.,|Bentz_etal 



2009ah ." closer to naive expectations from photoionization. The intrinsic scatt er of the R — L 



relation is estimated to be ~ 0.15 dex (~ 0.1 1 dex with the best quality RM data, Peterson 201 0|). 
The latest version of the R - L relation based on Hy6 RM measurements is ( Bentz et al.. 2009ai) : 



l°gr"^71 — = -21.3 +0.519 log ■ — . (3) 

light days ergs^' 



The tightness of the R - L relation has led to suggestions to use this re lation as an absolute 



luminosity indicator to use quasars as a cosmology probe (e.g.. lWatson et al. 2011; Czernv et al 
although RM measur ements of BLR siz es and quantification of the R - L relation beyond 
z ~ 0.3 are yet to come (e.g., Kaspi et ani2007 ). 



Single-epoch (SE) virial BH mass estimators The observed R - L relation provides a much 
less expensive way to estimate the size of the BLR based on the luminosity of the quasar. Sub- 
sequently this relation has been used to develop the so-called "single-epoch virial black hole mass 
estimators'0 (SE virial masses or SE masses in short hereafter): one estimates the BLR size from 



^It is also important to recognize that AGN SEDs can vaiy significantly from object to object, and some of these 
variances in SED must introduce certain scatter in the observed R - L relation based on measurable continuum instead of 
the ionizing continuum, and may cause systematic changes of BLR structure with SED properties (which might affect the 
virial coefficient /). 

'Such methods are also known as BH mass scalin g methods (e.g., Vesteraaard et al. 201 1). Occasionally, this method 
is referred to as the "photoionization method" (e.g.. |Petersonll201lLSalviander & Shields,, 20121) . This may be a little 
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the measured quasar luminosity using the R - L relation, and the width of the broad emission 
line, which are then combined to give an estimate of the BH mass using calibration coefficients 
determined from the sample of AGNs with RM mass estimates. Specifically, these estimators 
take the form: 

logl — I = a + felog|— rr^^ r| + clog( ^ , ) , (4) 

XMqI ^\1044ergs-i/ ^Vkms-'/ 

where L and W are the quasar continuum (or line) luminosity and width for the specific line, 
and coefficients a, b and c are calibrated against RM AGNs. In Eqn. Q, the coefficient on 
line width is usually taken to be c = 2, as expected from viral motion. Other values of c are 
suggested, however, depending on the definition of line width, and will be discussed further in 
|3] Based on the general similarity of quasar SEDs (both continuum and line strength), dif- 
ferent luminosities have been used, including continuum luminosities in X-ray, restframe UV 
and optica l, as well as line luminosities, in various versions of these single-epoch virial estim- 
ators (e.g.,'Vestersaard"2002^, 'McLure & Jarvis"2002l "iVIcLure & DunloD''2004'; 'Wu et al."2004'; 
Greene & Ho 2005; Vesteraaard & Peterson 2006; Kollmeier et al. 2006; Onken & Kollmeier 2008; 
Wang e t al. 2009; Vestergaard & Osmer2009; Greene et al. 2010b(; IRafiee & Halll201 IbtlShen et al] 
201lt IShen & Liu,2012 ; Trakhtenbrot & Netzer ,2012) . In general continuum luminosities 



are 



preferred over line luminosities given their tighter correlations with BLR size, but in some cases 
line luminosit ies are preferred wh ere the continuum may be significantly contaminated by host 
starli ght (e.g., Greene & Hol2005 ), or by the nonthermal emission from a jet in radio-loud objects 



(e.g.. lWuetal.ll2004 . As for the choice of line width, both FWHM and line dispersion (o-\me) are 



utilized in these calibrations. 

The uncertainty of these various single-epoch virial estimators can be inferred from the re- 
sidual s in the calibrations against the RM masses, and is esti mated to be on the order of ~ 0.5 dex 
(e.g., McLure & Jarvisll2002 ; Vestergaard & Petersonll200"6 ). This is similar to the uncertainty of 



RM masses, and can hardly be smaller since this method is rooted in the RM technique. Reasons 
and consequences of such substantial mass uncertainties will be elaborated in ^ From now on I 
will refer to RM and SE masses collectively as virial BH masses. 

The virial estimators (RM and SE) currently are the best method in estimating quasar BH 
masses. Therefore after a brief discussion on alternative methods 02. 2b . I will focus on these 
virial estimators in the rest of the review. 



2.2 Other Methods to Estimate Quasar BH Masses 



There are several other methods to estimate the mass of quasars. They are much less popular than 
the RM method and its extension, the SE virial method. Nevertheless there is certain merit in 



ambiguous, since in practice tliis empirical method is based on RM results rather than photonionization calculations, even 
though the observed R - L relation is consistent with naive photoionization predictions. The "photoionization method" 
better refers to those that estimate the BLR size using photoionization arguments (see j|2.21 . 
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further developing some of these methods, for instance, to provide complementary mass estimates 
and consistency checks. Therefore I give a brief discussion on these alternative methods. 

Photoionization method Historically Dibai was the first to sys tematically measure BH masses 
for quasar samples since the 1970s (e.g.. lDibaJl977[ll980l l 19841) . Adopting Woltjer's po stulation 



that th e BLR gas is in virial equilibrium in the gravitational potential of the central BH dWoltier 



19591), Dibai used Eqn. (|2| to estimate the BH mass, where the width of the line is used to indicate 
the virial velocity. To estimate the BLR size R, Dibai used the photoionization argument: 



Um = —R^j{n,T,)ev , 



(5) 



where L(Hfi) is the luminosity of the li/3 line, jirtgTe) is the volume emissivity in the H/3 line from 
photoionized gas, and ey is the volume filling factor of BLR clouds. Adopting constant values of 
He ~ lO' cm"-*, Ti; ^ lO'* K and ey ~ lO"-', Dibai estimated BH masses for more than ~ 70 nearby 
Seyfert 1 galaxies and quasars, and made the first plot of the distribution of AGNs in the mass- 
luminosity plane. Technically speaking, Dibai's method is also a single-epoch method, and it has 
an effective/?-/, relation ofR oc L(HjS)^^^, shallower than the observed /?-L relation. Despite the 
simplifications in Dibai's approach (sometimes unphysical), many of his BH mass estimates of 
local AGNs are consistent with today's RM masses to within 0.3 dex (e.g., Bochkarev & Gaskell 
I2OO9I) . Dibai's method a lso motivated some late r quasar BH mass estimations based on the same 



argument. For example. IWandel & Yahil (11985 ') used the same method with modifications to the 
volume filling factor, and derived a radius-luminosity relation R oc L(Hy6)'^^. 

Along a completely independent path, photoionization argume nts based on the ion ization 
parameter U were used to estimate the BLR size (e.g., Netzer 1990; Wandel et alj|l999|) . build- 



Davidson 



ing on earlier development of photoionization equilibrium theory in the 1970s (e.g 
I972I : iDavidson & Netzenll979i) . This approach emphasizes more on the role of the ionization 
continuum, and provides an intuitive understanding of the observed R - L relation. But just as 
Dibai's method, all these photoionization-based methods require assumptions or indirect con- 
straints on the physical conditions of the BLR gas (such as density, covering factor, etc) to infer 
the BLR size, therefore today they are not as popular as the more empirical, but more accurate 
RM-based methods discussed in 32.11 

Accretion disk model fitting (SED fitting) Another method to infer the mass of the BH is by 
fitting the SED of quasars. T he development of accretion d isk theory over the last four decades 
(for a recent review, see, e.g., lAbramowicz & Fragile l l2013l) h as enabled predictions of the emit- 
ting continuum spectrum of accreting BHs. By fitting the observed quasar continuum SED, one 
can constrain the model parameters (such as BH mass, accretion rate, BH spin, inclination) with 
adequate accretion di sk models. Many studies have used this SED fitting method to infer BH 
masses in AGNs (e . g..lMalkanll983l:[Sun & Malkanll989l:IWandel & Petrosianll988H Laodll990t 
Rokaki et alJll992l: iTripp et al.lll994 iGhiselHni et al."2010'; 'Caldero ne et al.ll2012h . usuallv as- 
suming a standard thin accretion disk model CShakura & Sunvaev..l973h . One main concern here 
is that standard accretion disk models, while can successfully produce broad-band features (such 
as the "big blue bump", e.g., [S hields 1978), do not yet have the capability to fully explain the 
AGN SED (e.g.. lKoratkar & Blae&.1999:.Lawrencel2012h . and the resulting BH mass constraints 
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may be sensitive to these deviations from standard accretion disk models. Given that there are 
parameter degeneracies and model assumptions/simplifications in the SED fitting procedure and 
the requirement for good multiwavelength coverage in UV-optical (where most of the disk emis- 
sion comes from), this method generally canno t prov ide an accuracy of better tha n a factor of ~ 
five in BH mass estimates at the moment (e.g.. lLaoii ll990: Calderone et al.ll2012h . But it would 
be interesting to compare this method with the virial methods 02.1b for larger samples. 

Microlensing in gravitationally lensed quasars - alternative routes to BLR sizes Since resolving 
the BLRs requires //as to tens of yuas angular resolution, RM will remain the primary method to 
measure the size of BLR in the next decade or two. Another indirect method to measure BLR size 
is via microlensing in gravitationally lensed quasars. The Einstein radius in the source (quasar) 
plane of a point-mass lens with mass M is r£ = -sjAGMDsDis /(Dic^), where Dj, Di and Dis 
are the angular diameter distances to the source, to the lens, and from the lens to the sourse. For 
typical values zl 
scale 9e 



0.5 and zs 



2 we have r£ a; 0.01 y/MjM^pc, corresponding to an angular 
fs/Ds ~ -\/M/Mq fj.as. This scale is comparable to the size of the BLR and accre- 
tion disk. Due to the relative transverse motion between the lens and source, large magnification 
can happen on short time scales (fcross - ysource/t^ j.) during caustic-cros sing in the source plane 



(for a primer on gravitational lensing, see, e.g., ISchneider et al. I ll992i) . Since the continuum 
emission (from the inner accretion disk) and BLR emission have diff'erent spatial scalefl mi- 
crolensing will cause differential variability for the continuum and different broad lines, which 
can be used to constrai n the size (and geometry) of the em i tting regions, such as the accretion 
disk and the BLR (e.g.Jlrwin et al.|ll989l: iLewis et al.l'l998': 'Popovi c et al.ll200ll: iRichards et al " 



20p4a; Morgan et al. 2010; Dai et al. 2010; Mosquera & Kochanek 2 01 ll:ISluse et alJ201 llboTa 



Guerras et al.ll2013h . In particular, the latest study by Guerras e t al.l (l2013l) found a R — L rela- 
tion based on microlensing BLR sizes, which is in reasonably good agreement with that based on 
RM. This technique is primarily applied to quasars that are already strongly lensed (with multiple 
images) for which the microlensing probability by stars/compact objects in the foreground lens 
galaxy is high, and the time delays between different images are known. The latter is important to 
rule out variability due to intrinsic AGN variability, which will often complicate the microlensing 
interpretation. 



Direct dynamical BH masses Although observationally challenging (given the overwhelming 
AGN continuum that dilutes the stellar absorption features and nongravitational forces on gas 
dynamics in the nucleus), there have been several attempts to get direct dynamical measurements 



of BH masses in Type 1 AGNs, using spatially resolved stellar kinem atics (e.g., iDavies et al 



2006; Onken et al. 2007 ) or gas kinematics (e.g.. Hicks & MalkanlE0 08) down to the sphere of 
influence, Rsi = GMsn/o'l, of the BH. The number of AGNs with reliable dynamical BH mass 
measurements is still small, and it would be important to obtain more dynamical mass measure- 
ments of AGNs to provide critical consistency checks on virial BH masses. 



The BLR and the accretion disk are probably not disjointed, i.e., some portion of the BLR gas could originate from 
the outer part of the accretion disk or from a wind launch ed from the disk, and high-ionization Unes could have differen t 
origins from low-ionization lines, as in some models (e.g.. lCollin-Souflfrin et al.ll988l:lElvisl2000l ; lRisaliti & Elvisl2010l) . 
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Other indirect methods These methods use correlations between BH mass and other measurable 
quantities, calibrated with known BH masses, to infer the BH mass in quasars. By virtue they can 
be applied to non-broad-line quasars as well, if the required quantity is measurable, but caution 
should be paid to possible systematics of each method. These correlations include the well-known 
BH-host scaling relations found for local inactive galaxies. For instance, the Mbh - cr» relation 
in local inactive galaxies is used to predict BH masses in local RM AGNs with cr, measurements 
(see ^3.1.21 ). Another example is the observed anti-correlations between BH mass (estimated 
with other methods) and X-ray variability properties of AGNs, such as the break frequ ency in the 
powe r spectral densities ( PSDs) of X-ray hght curves or variabi lity amplitude (e.g.. IPapadakisI 



2004 lO'Neil etalJl2005l: iMcHardv etaDbOOd IZhou et al.ll201Q) . which have the potential to 
provide independent BH mass estimates to within a factor of a few to the RM masses. The 
reliability and systematics of the X-ray variability method, however, are yet to be explored with 
larger AGN samples with known BH masses. 



3. Caveats, Uncertainties, and Biases of Virial BH Masses 



3.1 Physical Concerns 



3.1.1 The virial assumption 



There are eviden ce support ing the virial assumption in RM in at least several AGNs (e.g..lPeter son & Wandel 
[1999, 2000; Onk en & Pete rson 2002; Kollatschnv 1 2003h . For these objects RM lags have been 
successfully measured for multiple lines with different ionization potentials (such as Hy6, Civ, 
Hen) and line widths, which are supposed to arise at different distances, as in a stratified BLR 
for different lines. The measured lags and line widths of these different lines fall close to the 
expected virial relation W oc R^'^l^, although such a velocity-radius scaling does not necessarily 
rule out other BLR models where the dynamics is not dominated by the gravity of the central BH 
(e.g., see discussions in Krolik 2001 ). A more convincing argument is based on velocity-resolved 
RM, where certain dynamical models (such as outflows) can be ruled out based on the d ifference 
(or lack thereof) in the lags from the blue and red parts of the line (e.g.. iGaskellll 19881) . On the 
other hand, non-virial motions (such as infall and/or outflows) may indeed be pr esent in some 
BLRs, as inferred from recent velocity-res olved RM in a handful of AGNs (e.g.. 



Dennev et al 



2009aHBentz etalJl2010l: iGrier et al.ll2013h . Fortunately, even if the BLR is in a non-virial state. 



one might still expect that the velocity of the BLR clouds (as measured through the line width) 
does not deviate much from the virial velocity. Thus using Eqn. ^ does not introduce a large 
bias, and in principle this detail is accounted for by the virial coefficient / in individual sources. 



A further test of the virial assumption on the single-epoch virial estimators is to see if the 
line width varies in accordance to the changes in luminosity for the same object. The picture 
here is that when luminosity increases (decreases) the BLR expands (shrinks), and the Une width 
should decrease (increase), given enough response time. This test is important, because if the 
Une width does not change accordingly to luminosity changes, the SB mass will change for the 
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Figure 1. A test of the virial assumption using two-epoch spectroscopy from SDSS for Hyfl (upper), Mg ii 
(middle) and CiV (bottom). Plotted here are the changes in line width as a function of changes in continuum 
luminosity (Lsioo, ^3000, and L1350 for Hy6, Mgn, and CiV, respectively) between the two epochs. The left 
column is for FWHM and the right column is for crii„e. The dots are for all objects with measurement S/N> 3 
at both epochs for both L and W (not for AlogL and A log W). These objects tend to cluster around zero 
values because the typical continuum and line variabilities of SDSS quasars are small. The red triangles are 
the median values in each A log L bin, where the eiTor bars indicate the uncertainty in the mean. A perfect 
virial relation would imply A log W = -0.25A log L, as indicated by the solid line in each panel. Note that 
I have neglected the chromatic nature of quasar variability, which would predict an even steeper relation 
between A log W and AlogL (see g3.1.1l for details). The low-redshift z g 0.7 SDSS quasars with median 
luminosity (log(L5ioo/erg s"')) = 10'*'"' show the expected virial relation between AlogL and AlogW, 
which is not the case for the high-luminosity SDSS quasars at z > 0.7 based on Mg 11 or CiV. 
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same object, introducing a luminosity-dependent bias in the mass estimates (see ^3.3.21 ). This 
test is challenging in practice, given the limited dynamic range in continuum variations and the 
presence of measurement errors. Nevertheless, in several AGNs with high-quality RM data, such 
anti- correlated va riations of li ne width and BLR size (or continuum luminosity) have been seen 
(e.g., Peterson et al. 2004; Park et all2012b ). once the lag between continuum and line variations 



is taken into account. While this lends some further support for RM and SE virial estimators, it 
should be noted that: 1) not all RM AGNs show this expected behavior, given insufficient data 
quality; 2) it makes a difference which line width measurements (i.e., FWHM vs cr, rms vs mean 
spectra) and which BLR size estimates (i.e., t vs continuum luminosity) are used. 

It is also not clear if the above results based on a few RM AGNs apply to the general quasar 
population. Fig.[T]shows a test of the co-variation of line width and continuum luminosity using 
thousands of SDSS quasars with spectra at two epochs (She, Shen, et al., in prep). While for the 
majority of these quasars the two epochs do not span a large dynamic range in luminosity, the 
large number of objects provide good statistical constraints on the average trend. In Fig. [T]the 
black dots are measurements for individual objects, and they cluster near the center because most 
quasars do not vary much between the two epochs. The measurement uncertainties on AlogL 
and A log W are large, so I bin the results in A log L bins and plot the medians and uncertainties 
in the median in each bin in red triangles. The measurement uncertainties in A log L and A log W 
are comparable for all three lines, but only for the low-luminosity and low-z (z < 0.7) H/? sample 
is the median relation consistent with the virial relation (the solid lines in Fig.[T]i. For the other 
samples at z > 0.7 based on Mg ii and Civ, the line width does not seem to respond to luminosity 
changes as expected from the virial relation. This difference could be a luminosity effect, but 
more detailed analyses are needed (She, Shen, et al., in prep). 

Another important point to make is that there is a well know n fact that quasar spectra get 
harder (bluer) as they get brighter (e.g., IVanden Berk et al. lllooi and references therein). This 



means that the variability amplitude in the ionizing continuum should be larger than that at longer 
wavelengths (i.e., the observed continuum). Thus we should see a somewhat steeper slope in the 
line width change versus the (observed) continuum luminosity change plot for a single object 
(e.g., iPeterson et al ] l2002h . This, however, would be in an even larger disagreement with the 



trends we see in Fig.[T] 



3.L2 The virial coefficient f 



To relate the observed broad line width to the underlying virial velocity (e.g., Eqn.|2l) requires the 
knowledge of the (emissivity weighted) geometry and kinematics of the BLR. In principle RM 
can provide such information, and determine the value of / from first principles. Unfortunately 
the current RM data are still not good enough for such purposes in general, although in a few 
cases alternative approaches have been invented lately to ac count for the effect of / in directly 



modeUng the RM data using dynamical BLR models (e.g.. iBrewer et al.ll20I It iPancoast et al 



2012h . Early studies made assumptions about the geornetry and stru cture of the BLR in deriv 



ing RM masses (e.g.. lNetzei]|l990t IWandel et al.lll999l; iKaspi et al.l 12000.) or SE virial masses 
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(e.g., McLure & DunlopI 2004 ). Now the average value of / is mostly determined empirically 
by requiring that the RM masses are consistent with those predicted from the Mbh - cr» relation 
of local inactive galaxies. Such an exercise was first done bv lOnken etalj (|2004 . who used 16 
local AGNs with both RM measurements and stellar velocity dispersion measurements to derive 
(/) a! 1.4 if FWHM is used, or (/) » 5.5 if cr\i„f. is used. Later this was repeated with new RM 
data (e.g.. iWoo et al.ll201Q) . who derived a similar value of (/) a: 5.2 (using criine)- 

However, in recent years it has become evident that the scaling relations between BH mass 
and bulge properties are not as simple as we thought: it appears that different types of galax- 
ies follow somewhat diffe rent scaling r elations, and th e scatter s eems to increase towards less 
massive systems (e.g..lHur2 008; Greene et al.l2008l:l&ah am 2008l: lGraham & Lil2009trHul2009 



Giiltekin et al.l2009[ Greene et al l2010bi:lMcConnell & Mal20li and references therein). There 
fore, depending on the choice of the specific form of the Mbh - cr* relation used and the types 
of galaxies hosting R M AGNs in th e cahb ration, the derived average / value could vary signi- 
ficantly. F or instance, Graham et al. (2011) de rived a { /) value that is is o nly half of the values 
derived bv lOnken eFal.. (.2004.) and ,Woo et al.1 (l2010l) . (Park et all (l2012bl) performed a detailed 
investigation on the effects of different regression methods and sample selection in determining 
the Mbh - cr* relation and in turn the (/) value, and concluded that the latter is the primary cause 
for the discrepancy in the reported (/) values. Given the small sample sizes of RM AGNs with 
host property measurements and the uncertainties in the BH-host scaling relations in inactive 
galaxies, the uncertainty of (/) is still ~ a factor of 2 or more, and will remain one of the main 
obstacles to estimate accurate RM (or SE) BH masses in terms of the overall normalization. One 
may also expect that the actual / value is different in individual sources, either from the diversity 
in BLR structure or from orientation effects (since the line width only reflects the line-of-sight 
velocity, see ^3.1.6l l. Thus using a constant / value in these RM masses and SE virial estimators 
introduces additional scatter in these mass estimates. 



Perhaps a more serious concern is the assumption that the BH-host scaling relations are the 
same in active and inactive galaxies. While there is a cle ar correlation between bulge properties 
and the RM masses in RM AGNs (e.g., Bentz et alj2009 cV it could be offset from that for inactive 
galaxies if the actual (/) value is different. Such a scenario is plausible if the BH growth and host 
bulge formation are not always synchronized. The only way to tackle this problem is to infer 
/ from directly constrained BLR geometry /kinematics with exquisite velocity-resolved RM data 
that map the line response (transfer function) in detail, and this must be done for a large number 
of AGNs to explore its diversity. 



3. 1 .3 FWHM versus line dispersion 



Both FWHM and triine are commonly used in SE virial mass estimates as the proxy for the virial 
velocity (when combined with the virial coefficient /). Both definitions have advantages and 
disadvantages. FWHM is a quantity that is easier to measure, less susceptible to noise in the 
wings and treatments of line blending than cr\^n^, while crjine is less sensitive to the treatment of 
narrow line removal and peculiar line profiles. Overall FWHM is preferred over cr\in^ in terms of 
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easiness of the measurement and repeat ability. As a -wne, measurements depend sensiti vely on data 



quality and different methods used (e.g. JPenney et al. 2009b; Rafiee & Hall 2011a b: lAssef et al 
1201 lb . the SE virial masses (e.g., Eqn. |4]i based on criine could differ significantly for the same 
objects. 

Physically one may ar gue crijne is more trus tworthy to use than FWHM, although the evidence 



to date is only suggestive. ICollin et al.l (120061) compared the virial products based on both crjine 
and FWHM with those expected from the Mbh - cr« relation, for 14 RM AGNs. All their line 
width measurements were based on the rms or mean spectra of the RM AGNs. They found that 
the average scale factor (i.e., the virial coefficient /) between virial products to the Mbh - o"» 
masses depends on the shape of the line if FWHM is used, while it is more or less constant if a-\iae 
is used. Based on this, they argue that crii„e is a better surrogate to use in estimating RM masses. 
Additionally, crwne measured in r ms spectra seems to fo llow the expected virial relation better than 
FWHM in some RM AGNs (e.g., Peterson et alj|2004 ). although such evidence is circumstantial. 



It is impo rtant to note that for a given line, the ratio of FWHM to crune is not necessarily a 
constant (e.g.. lCollin et al. 1 12004 |Petersonll201 iL but cf., Decarli et al. 2008a), while a Gaussian 
line profile leads to FWHM/crii„e ^ 2.35. For H/?, FWHM/criine seems to increase when the 
line width increases. This might be related to th e Populations A and B sequences developed by 
Sulentic and collaborators (ISulentic et al .l2000ah . which is an extension of earlier work on the cor- 
relation space of AGNs (the so-called "eigenvector 1", e.g., Bor oson & Greenlll992tlWang et al 
[1996). A direct consequence is that there will be systematic differences in Mse whether FWHM 
or cTiine is used for the same set of quasars, especially for objects with extreme l ine widths. In 
general a "tilt" between the FWHM and crune-based virial masses is expected (e.g.. lRafiee & Hall 
I2OI la.b ). Currently directly measuring criine from single-epoch spectra is much more ambiguous 
and methodology-dependent than measuring FWHM. If one accepts that crune is a more robust 
virial velocity indicator, it is possible to c onvert the measured FWHM to crijne using the relation 
found for high S/N data (e.g..l Collin et al..,20 06). or empirically determine the depen dence of SE 



mass on FWHM (i.e., coefficient c in Eqn.|4]l using RM masses as calibrators (e.g.. lWang et al 
20091) . which generally leads to values of c < 2. 



The choice of l ine width indicator s is still an open issue. 

(I2OO6"), 



arguments in, e.g., Collin et al. 



It will be important to revisit the 
using not only more but also better-quality RM data, 
as well as to investigate the behaviors of FWHM and crune (and perhaps alternative line width 
measures) for large quasar samples. 



3 . 1 .4 Broad line profiles 

As briefly mentioned in ^2.11 part of the reason that we are struggling with / and line width 
definitions is because of the simplifications of a single BLR size and using only one line profile 
characteristic to infer the underlying BLR velocity structure. If we have a decent understanding 
of the BLR dynamics and structure (geometry, kinematics, emissivity, ionization, etc.), then in 
principle we can solve the inverse problem of inferring the virial velocity from the broad line 
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profile. Unfortunately, the detailed BLR properties are yet to be probed with velocity-resolved 
reverberation maps, and the solution of this inverse problem may not be unique (e.g., different 
BLR dynamics and structure may produce similar line profiles). 



Nevertheless, there have been efforts to model the observed bro ad line profiles with simple 

BLR models. The best know n example is the disk-emitter model (e.g., Chen et al.l 19891 : Eracleous & HalpenJ 
1994 ; Eracleous et al. 1995h . where a Keplerian disk with a turbulent broadening component is 
used to model the double-peaked broad line profile seen in ~ 10 - 15% radio-loud quasars (and 
several percent of radio-quite quasars). The line profile then can place constraints on certain 
geometrical parameters, such as the inclinati on of the disk, thus has relevance in the / value for 
individual objects (e.g.. iLa Mura et al. I l2009t) . Another example is using simple kinematic BLR 
mod els to explain the trend of the line shape parameter FWHM/criine as a function of line width 
re.g..l Kollatschnv & Zetz]|[2oT l M2013h . as mentioned earlier in 33.1.31 These authors found that 
a turbulent component broadened by a rotation component can explain the observed trend of 
line shape parameter, and their model provides conversions between the observed line width and 
the underlying v irial (rotational) velocity. More complicated BLR models can be built (e.g.. 
Goad et al. 20121) . which has the potential to underpin a physical connection between the BLR 
structure and the observed broad line characteristics. While all these exercises are worth further 
investigations, it is important to build self-consistent models that are also verified with velocity- 
resolved RM. 



3 . 1 .5 Effects of host starlight and dust reddening 



The luminosity that enters the R - L relation and the SB mass estimators (Eqn. |4) refers to the 
AGN luminosity. At low AGN luminosities, the contamination from host starlight to the 5100A 
luminosity can b e significant. This m otivated the alternative uses of Balmer line luminosities 
in Eqn. Q (e.g., Greene & Hoi 2005b . Using line luminosity is also preferred for radio-loud 
objects whe re the continuum may be severely c onta minated by the nont hermal emission from 
the jet (e.g.. lWuet"ai]|2004 . iBentz et alj ( l2006h and lBentz et all ( l2009al) showed that properly 
accounting for the host starlight contamination at optical luminosities in RM AGNs leads to a 
slope in the R - L relation that is closer to the naive expectation from photoionization. Similarly, 
using host-c orrected L5100 can le ad to reduced scatter in the HjS-Lsioo SE calibration against RM 
AGNs (e.g.. lShen & Kellvll2012h . 



The average contribution of host starlight to L5100 has been quantified bv lShen et alJ (1201 Ih . 
using low-redshift SDSS quasars. They found that significant host contamination (> 20%) is 
present for log Lsioo.totai < lO"*'*^ erg s"', and provided an empirical correction for this average 
contamination. Variations in host contribution could be substantial for individual objects though. 



For UV luminosities (Lsooo, ^^1350 or L1450), the host contamination is usually negligible, al- 
though may be significant for rare objects with excessive ongoing star formation. A more serious 
concern, however, is that some quasars may be heavily reddened by dust in ternal or external to 
the host. The so-called "dust-reddened" quasars (e.g.. iGhkman et alj|2007l) have UV luminos- 
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ities significantly dust attenuated, and corrections are required to measure their intrinsic AGN 
luminosities. It is possible that optical quasar surveys (such as SDSS) are missing a significant 
population of dust-reddened quasars. 



3. 1 .6 Effects of orientation and radiation pressure 



If the BLR velocity distribution is not isotropic, orientation effects may affect the RM and SE 
mass estimates. Specific BLR geometry and kinematics, such as a flattened BLR where the orbits 
are confined to low latitudes, will lead to orientation-dependent line width. Some studies report a 
correlatio n between the broad line FWHM and the source orientation inferred from radio proper- 
tied! (e.g., IwmT&Brownilllii ; lJarvis_ & McLure 2006)^in favor of a flattened BLR geometry. 
Similar conclusions were achieved in Decar li et al.l (l2008al) based on somewhat different argu- 
ments. Since we use the average virial coefficient (/) in our RM and SE mass estimates, the true 
BH masses in individu al sources may be over- or underestimated depending on the actual incli n- 
ation of the BLR (e.g.. lKrolikll200ll:lDecarU et alj|2008al: iFine et al.ll20Tl1; feunnoe et alj|2013lfl 



The distributions of broad line widths in bright quasars are ty pically log- normal, wit h dispersions 
of ~ 0.1 - 0.2 dex over ~ 5 magnitudes in luminosity (e.g.. lShen et al ]|2b08 a; Fine et al.ll2008l 



20I0h . A thin disk-like BLR geometry with a large range of inclination angles cannot account 
for such narrow distributions of line width, indicating either the inclination angle is limited to a 
narrow range for Type 1 objects, and/or there is a significant random velocity component (such as 
turbulent motion) of the BLR. This Umits the scatter in BH mass estimates caused by orientation 
eff'ects to be < 0.2 - 0.4 dex. 



So far we have assumed that the dynamics of the BLR is dominated by the gravity of the 
central BH. The possible effects of radiation pressure, which al s o has a oc R^^ dilution as gravity, 
on the BLR dynamics have been emphasized by, e.g., iKrolikI ( 120011) . On average the possible 
radiation effects are eliminated in the empirical calibration of the (/) value (see ^3.L2l l. but 
neglecting such effects may introduc e scatte r in individual sources and luminosity-dependent 
trends. Most recently iMarconi et al. I (l2008h modified the virial mass estimation by adding a 
luminosity term: 

MBH.M0B=/^+^(y^54^)^O, (6) 

where the last term describes the effect of radiation pressure on the BLR dynamics with a free 



^Some recent studies (e.g. jFine et al.ll2bl iHRunnoe et alj|2013i) argue that the dependence of hne FWHM on source 
orientation is different for low-ionization and liigh-ionization lines, such that the CIV-emitting gas velocity field may be 
more isotropic than H/3 and Mg ii. 

*Of some relevance here is the interpretation of the apparently small BH masses in a sub-class of Type 1 AGNs 
called nan'ow-line Seyfeit Is (NLSls), where the H/J FWHM is narrower than 2000 km s~* alon g with other unusual 
properties (such as strong iron emission and weak [O in] emission). Some argue (e.g..l Decarli et a l. 2008b) that NLSls 
are preferentially seen close to face-on, hence their virial BH masses based on FWHM are underestimations of true 
masses. However, NLSls also differ from normal Type 1 objects in ways that are difficult to explain with orientation 
effects (such as weak [Oiii] and strong X-ray variability). Orientation may play some role in the interpretation of NLSls 
(especially for a minority of radio-loud NLSls), but is unlikely to be a major factor. 
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parameter g. By allowing this extra term, 'Marc oni et al.l (120081) re-calibrated the RM masses 
using the Mbh - o"» relation, and the SE mass estimator using the new RM masses. This approach 
improves the rms scatter between single-epoch masses and RM masses, from ~ 0.4 dex to ~ 0.2 
dex, and rem oves the slight systematic trend of the SE mass scatter with RM masses seen in 



Vestergaard & Peterson (2006). However, it is also possible that the reduction of scatter between 
the SE and RM masses is caused by the addition of fitting freedoms. Since the intrinsic errors on 
the RM masses are unlikely to be < 0.3 dex, optimizing the SE masses relative to RM masses to 
smaller scatter may lead to blown-up errors when appl y the optimized scali ng relation to other 
objects. It would be interesting to split the RM sample in lMarconi et al. I (l2008h in half and use one 
half for calibration and the other half for prediction, and see if similar scatter can be achieved in 
both subsets. The relevance of radiation pressure is also questioned by Netzer (2009), who used 
large samples of Type 1 and Type 2 AGNs from the SDSS to show that the radiation-pressure 
corrected viral masses lead to inconsistent Eddington ratio distributions in Type Is and Type 
2s, even though the [O m] luminosity distribution is consistent in the two samples. However, 



ougnme t 

Marconi et al. [ (120091) argues that the diff'erence in the "observed" Eddington ratio distributions 
does not mean that radiation pressure is not important, rather it could result from a broad range of 
column densities which are not properly described by single values of parameters in the radiation- 
pressure-corrected mas s formula. These stu dies then revealed that using the simple corrected 
formula as provided in iMarconi et al.l (120081) does not provide a satisfactory recipe to account 
for radiation pressure in RM or SE mass estimates, and the relevance of radiation pressure and 
a practical method to co rrect for its effect are therefore still under active investigations (e.g., 
Netzer & Marzianill2010l) . 



3. 1 .7 Comparison among different line estimators 



There are both low-ionization and high-ionization broad lines in the restframe UV to near-infrared 
of the quasar spectrum. Despite different ionization potential and probably different BLR struc- 
ture, several of them have been adopted as SE virial mass estimators. The most frequently used 
line-luminosity pairs include strong Balmer lines (Ha and H/3) with Lsioo or ijjQ, jjyj, Mg ii with 
L3000, and Civ with Logo or Lhsq. Hydrogen Paschen lines in the near-IR can also be used if 
such near-IR spectroscopy exists. 



There have been SE calibrations upon specific lines against RM masses, or against SE masses 
based on another line. Comparisons between different SE line estimators using various quasar 
samples are often made in the literature: s ome claim consistency, while others report discrep- 
ancy. As emphasized in lShen et al.l (l2008al) . it is important to use a consistent method in meas- 
uring luminosity and line width with that used for the calibrations if one wants to make a fair 
comparison using externa l samples. Failure to do so may lead to unreliable conclusions (e.g.. 



Dietrich & Hamannll2004 . 



The continuum luminosities at different wavelengths and several line luminosities are all cor- 
related with each other, with different levels of scatter Fig. |2] shows some correlations between 
different continuum luminosities using the spectral measurements of SDSS quasars from .Shen et al. 
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teOllh . To compare L1350 and L5100 directly, one needs either UV+optical or optical+near-IR to 
cover both restframe wavelengths. Fig. |3](left) shows such a coniparis^njTom a recent sample 
of quasars with optical spectra from SDSS and near-IR spectra from iShen & Liul (l2012i) . which 
probes a higher luminosity range Lgioo > lO'^^^ergs"' than the SDSS sample. Correlations 
between these luminosities are still seen at the high-luminosity end. For the SDSS quasar popu- 
lation, different luminosities correlate with each other well, but this may be somewhat affected by 
the optical target selection of SDSS quasars that may preferentially miss dust-reddened quasars 
(see ^3.1.5b . In other words, the intrinsic dispersion in the UV-optical SED may be larger for the 
general quasar population. For instance. lAssefet al.lfcOl ll) found a much larger dispersion in the 
L1350/L5100 ratio for a gravitationally lensed quasar sample, which is selected differently from the 
SDSS. This large dispersion in the ^1350/^^5100 ratio will lead to more scatter between the H/? and 
Civ based SB masses. 



It is also important to compare the widths of different lines. Since Hy6 is the most studied line 
in reverberation mapping and the R — L re lation was measured using BLR radius for HyS (e.g.. 



Kaspi et alJ[2000ll2005l:lB"entz et alJl2009 aV it is reasonable to argue that the SB mass estimators 



based on the Balmer lines are the most reliable ones. The width of the broad Hot is well correlated 
with that of the broa d H/? and therefore it provides a good substitution in the absence of HyS (e.g., 

s of Mgii are found to correlate well with those of the Balmer 
"Salviander et al. 2007: McGill et a l. 2008; Shen et al. 2008a, 201 1; Wang et al. 2003; 



Greene & 



lines (e.g.. 



Vestergaard et al.ll201 It IShen & Liull2012[ see Fig. |2]for a comparison based on SDSS quasars). 



But such a correlation may not be linear: despite different methods to measure line widths, most 
recent studies favor a slope shallower than unity in the correlation between the two FWHMs (e.g., 
see Fig.|2ll. Given this correlation it is practical to use the Mg 11 width as a surrogate for HyS width 
in a Mg ii-based SB mass estimators, and some recent Mgn calibrations can be found in, e.g., 
IVester gaard & Osmer (2009): Shen & Liu (2012); Trakhtenbrot & Netzer (2012). However, one 
intriguing feature regarding the Mg 11 line is tha t the distribution of its line widths s eem to have 
small dispersions in large quasar samples (e.g., 'S hen et aL II2OO8 a; Fine e t al.ll2008l) . It appears 
as if the Mgn varies at a less extent compared with HjS (cf.,3Voo 2008, and references therein). 
It is also recently argued that for a small fraction of quasars (~ 10%) in the NLSl regime (e.g., 
small Hy6 FWHM and strong Feii emission), Mg 11 may have a blueshifted, non-virial component, 
and a n overall larger FWHM than HyS, that will bias the virial mass estimate (e.g., Marziani etjal] 



2013b . This i s consistent with the general trend found between Mg 11 and H j6 FWHMs using SDSS 



quasars (e.g., Wang et aDl2009 : Shen et aDl201 1 : Vestergaard et al.ir201 1 ). and may be connected 
to the disk wind scenario for Civ discussed below. 



The correlation between Hy8 (or Mg 11) a nd Civ widths is more controversial. While some 
clairn that these two do not correlate well (e.g.. iBachev et al l2004UBa skin & Laoi<2005"l: iNetzer et al. 
20071 IShen et al.l2008al: iFine et al.ll20 1ot lShen & Liul2012l:lT rakhten brot & Netzedl2012l). o thers 



claim there is a significant correlation (e.g. .Vestergaard & Petersonl2006t Assef et al .I2OII V Fig 
|3] (right) shows a compilation of Civ and H/3 FWHMs from the literature, which are derived for 
quasars i n different luminosities and red shift ranges. Only the low-luminosity (and low-z) RM 
sample in IVestergaard & Peterson (|2006^ shows a significant correlation. It is often argued that 
sufficient data quality is needed to secure the Civ FWHM measurements, although measurement 
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Figure 2. Comparisons between different continuum luminosities and line FWHMs, using SDSS quasar 
spectra that cover two lines. Shown here are the local point density contours. Measurements are from 
IShenetaLlJioT 3). The upper panels show the correlations between continuum luminosities, and the bottom 
panels show the correlations between line FWHMs. While the Mg n FWHM correlates with Hfi FWHM 
reasonably well, the correlation between the CiV FWHM and Mg ii FWHM is poor (also see, e.g. JShen et alJ 
l2008j : lFine et al.ll2008l.l2oTol) . 
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Figure 3. Left: correlations between different luminosities using the quasar sample in IShen&Liul(l2012h . 
which covers all four lines (CiV, Mgii, H/?, Ha) in the same object, for the high-luminosity regime 
isioo > lO'*' '* erg s"'. The solid lines are the bisector linear regression results using the BCES estimator 
("e.g.. lAkritas & Bershadvl[l996h . and the dashed lines indicate a linear correlation of unity slope. Right: 
comparison between CiV FWHM and Yip FWHM using different samples from the literature llShen & Liul 
j20ll 60 objects; SL12), lAssef et alj j201lL 9 obj ects; All), IVestergaard & PetersonI (l2006l 21 objects; 
VP06), iNetzer et alJ ( l2007l 15 objects; N07), and iDietrich et all ( l2009 l. 9 objects; D09)]. Only for the 
low-redshift and low-luminosity VP06 sample is there a significant correlation between the two FWHMs. 



errors are unlikely to account for all the scatter seen in the comparison between Civ and HjS 
FWHMs - the correlation between the two is still considerably poorer than that between Mgn 
and Hj6 FWHMs for the samples in Fig.[3]when restricted to high-quality data. IShen & Liul ( 120121) 
suggested that the reported strong correlation between Civ and H/? FWHMs is probably caused 
by the small sample statistics, or only valid for low-luminosity objects. 



The high-ionization Civ line also differs fr om low-ionization lin es such as Mgii and the 
Balmer lines in many ways (for a review, see ISulentic et al.ll2000bl) . Most notably it shows 



a prominent blueshift (typically hundreds, up to thousands of kms ') with re spect to the low 



ionization lines (e.g.. lGaskelllll982l: iTvtler & Fanlll992t [Richards et al.l 120021) . which becomes 



more prominent when luminosity increases. There is also a systemic trend (albeit with large 
scatter) of increasing Civ FWHM and li ne asymmetry when the Civ blueshift increases, a trend 
not present for low-ionization lines (e.g., She n et al ] l2008al l 120111) . The Civ blueshift is predom- 
inantly believed to be an indication of outflows in some form, and integrated in the disk-wind 
framework discussed below (but see lGaskellll2009l for a different interpretation). These prop- 
erties of Civ motivated th e idea that Civ is l ikely more affected by a non-virial component than 
low-ionization lines (e.g.. lShen et al. 2008a|). probablv from a radiatively-driven (and/or MHD- 
driven) accre tion disk wind (e.g., iKonigl & Kartidll994t iMurrav et al.lll995l: iProga et al. I l2000t 



Everettll2005h . especially for high- luminosity objects. A generic two-component rnodel for the 



Civ emission is then implied (e.g., ICollin-Souffrin et al.lll988t iRichards et al.ll201 It I Wang et al 
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2011). A similar argument is proposed by bennev" ( 2012 ) based on the Civ RM data of local 
AGNs, where she finds that there is a component of the Civ line profile that does not rever- 
berate, which is likely associated with the disk wind (although alternative interpretations exist). 
This may also explain the poorer correlation between Civ width and H/? (or Mg ii) width for 
more luminous quasars, where the wind component is stronger (see further discussion in §l3.1.9t . 
Therefore Civ is likely a biased virial mass estimator (e.g., Baskin & Laorll2005l : ISulentic et alJ 
20071: iNetzer et al]l2007l: IShen et al.ll2008al: iMarziani & Surenticii2012l and references therein). 



Although in principle certain properties of Civ (such as line shape parameters) can be used to 
infer the Civ blueshift and then correct for the Civ-based SE mass, such corrections are difficult 
in practice given the large scatter in these trends and typical spectral quality. Proponents on the 
usage of Civ line often emphasize the need for good-quality spectra and proper measurements of 
the line width. But the fact is Civ is indeed more problematic than the other lines, and there is 
no immediate way to improve the Civ estimator for high-redshift quasars, although some recent 
works are show ing some promising trends that may be used to improve the Civ estimator (e.g., 
DennevlEoil) . 



There have also been proposals for using the Cm], A i m, or Sim] l i nes in replacement of 



Civ (e.g., iGreene et al.1 1201 Oat IMarziani & Sulenticll2012l) . IShen & Liul (120 12h found that the 



FWHIVis of Ci v and Cm] are c orrelated with each other, and hence Cm] may not be a good line 
either (also see lHo et al.ll2012h . On the other hand. Aim and Sim] are more difficult to measure 
given their relative weakness compared to Civ and Cm] as well as their blend nature, hence are 
not practical for large samples of quasars. Another possible Une to use is Lya. Although Lya is 
more severely affected by absorption, intrinsically it may behave similarly as the Balmer lines. 
Such an investigation is ongoing. 



To summarize, currently the most reliable lines to use are the Balmer lines, although this 
conclusion is largely based on the fact that these are the most studied and best understood lines, 
and does not mean there is no problem with them. Mg ii can be used in the absence of the 
Balmer lines, although the lack of RM data for Mg ii poses some uneasiness in its usage as a 
SE estimator Civ has local RM data (though not enough to derive a R - L relation on its own), 
but the application of Civ to high-redshift and/or high-luminosity quasars should proceed with 
caution. In light of the potential problems with Civ, efforts have been underway to acquire 
near-IR spectroscopy to study the high-z quasar BH masses using Mg ii and Balmer lines (e.g., 
Shemmer et allbOO^ iNetzer et al.llioOTl IMarziani et al.ll2009HDietrich et al.ll2009tjGreene et al " 



2010al:lTrakhtenbrot et al. |201 lllAssef et al.l201 lllShen & Liul2012l:lHo et alj2012[|Matsuoka et al 
2013h 



3.1.8 Effects ofAGN variability on SE masses 

Quasars and AGNs vary on a wide range of timescales. It is variability that made reverberation 
mapping possible in the first place. One might be concerned that the SE masses may subject 
to changes due to quasar variability. Several studies have shown, using multi-epoch spectra of 
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quasars, that the scatter due to luminosity changes (and possibly corresponding change s in line 



width) does not introduce significant (> 0.1 dex) scatter to the SE masses (e.g., Wilh ite et al 
2007 : Dennev et al. 2009bt IPark et^l]|2012al) . This is expected, since the average luminosity 
varia bility a mplitude of quasars is only ~ 0.1 - 0.2 magnitude ove r month-to-year timescales 



(e.g., ISesaj et al. 2007; MacLeod et al. 2010l MacLeod et alJl2012l) . thus the difference in SE 



masses from multi-epochs will be dominated by measurement errors (in particular those on line 
widths). 



However it is legitimate to consider the consequence of uncorrected stochastic variations 
between line width and luminosity on SE masses, whether or not such uncorrelated variations are 
due to actual physical effects, or due to improper measurements of the continuum luminosity and 
line widths. Examples are already given in ^3.L1I and more detailed discussion will be provided 



3 . L9 Limitations of the RM AGN sample 



Last but not least, the current RM sample is by no means representative of the general quasar/AGN 
population. It is a highly heterogeneous sample that poorly samples the high-luminosity regime of 
quasars, and most objects are at z < 0.3. This alone calls into question the reliability of extrapola- 
tions of locally-caUbrated SE relations against these RM AGNs to high-z and/or high-luminosity 
quasars. 



The distribution of the RM AGNs in the spectral par ameter space of quasars is also highly 

biased relati ve to the general population . Richard s et al. (2011) developed (building on earlier 

ideas by, e.g.JComn-Sou ffrin et al.1 1 QSilMurrav et al..l995:.Proga et al.l2000[lElvisl2000l:lLeighlv & Moore 
|2004; Leighlv 2004) a generic picture of two-component BLR structure for Civ, composed of a 
virial component, and a non-virial wind component which is filtering the ionizing continuum from 
the inner accretion disk. This generic picture is able to explain, phenomenologically, many char- 
acteristics of the continuum and Civ line properties, such as the Civ blueshift and the Baldwin 
effect (i.e., the anti-correlation between Civ equivalent width and adjacent continuum luminosity, 
Baldwinll977 ). Fig.|4]shows the distribution of RM objects in the parameter space of Civ spectral 
properties, where most RM AGNs occupy the regime dominated by the virial component. Part 
of this is driven by luminosity, since more luminous quasars have on average larger Civ blueshift 
03. 1.7b . It will be important to ex plore this under-re presented regime with Civ RM at high- 
redshift, which has just begun (e.g. jKaspi et al.ll2007l) . Although this is an immediate concern 
for Civ, iRichards et al.l (1201 ll) made a fair argument that the BLR properties for Hj6 and Mg ii 
may also be biased in the RM sample relative to all quasars, if the non-virial wind component is 
also affecting the BLR of Hj6 and Mg ii by filtering the ionizing continuum. 



To date most of the RM lag measurements are for H ,6, and lag measurements are either 
lacking for Mg ii (but see lMetzroth et al.lE006l: IWodbOOSi and references therein, for Mg ii RM 
attempts and tentative results) or insufficient for Civ to derive a direct R - L relation based on 
these two UV lines. The total number of RM AGNs is also small, ~ 50, not enough to probe 
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CIV Blueshift [km/s] 

Figure 4. An updated version of Fig. 18 in lRichards et alj (1201 ih . showing the biased distribution of the 
local RM AGNs in the parameter space for CiV (blueshift relative to Mg ii versus the rest equivalent width). 
The contours and dots are 1.5 < z < 2.2 SDSS quasars, and the blue filled circles are RM AGNs. The black 
and red contours show the results for radio-quiet and radio- loud populations respectively. Most of the RM 
AGNs occupy quadrant I, where the CiV line is dominated by the virial component in the two-component 
picture in lRichards et al.l (1201 ll) . In this picture, quadrant IV is for CiV lines dominated by the non-virial 
wind component. The average quasar luminosity increases from quadrant I to quadrant IV. Figure courtesy 
of G. Richards. 
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the diversity in BLR structure and other general quasar properties. The current sample size and 
inhomogeneity of RM AGNs pose another major obstacle to develop precise BH mass estimators 
based on RM and its extension, SE virial methods. 



3.2 Practical Concerns 



3.2.1 How to measure the continuum luminosity and line widths 



Usually the continuum and line properties are measured either directly from the spectrum, or 
derived from fits to the spectrum with some functional forms for the continuum and for the 
lines. Arguably functional fits are better suited for spectroscopic samples with moderate to low 
spectral quality. As briefly mentioned earlier (§ 13.1.7) . it is essential to measure the continuum 
and line width properly when using the existing SE calibrations. Different methods sometimes 
do yield systematically different results, in particular f or the line width measure ments. Some 
studies fit the broad lines with a single component (e.g.. lMcLure & Dunlopll2004l) . while others 
use multiple components to fit the broad line. But if one wants to use the calibration in, e.g., 
McLure & Dunlopl ( 2004 ). then it is better to be consistent with their fitting method. Some com- 
parisons bet ween the broad lin e widths from different fitting recipes can be made using the catalog 
provided in lShen et al. ( 201 1 ). Take Hy6 for example, since this broad line is not always a single 
Gaussian or Lorentzian, the line widths from the single-component and multiple-component fit 
could differ significantly in some cases. 



The detailed descrip t ion of spectral fitting procedure can b e found in many papers (e. g . , 
McLure & DunlopI 120041: iGreene & Hoi 120051: IShen et al.ll2008a[ |201 ll: IShen & Liull2012h . In 
short, the spectrum is first fit with a power-law plus an iron emission templat^ in several spectral 
windows free of major broad lines. The best-fit "pseudo-continuum" is then subtracted from 
the original spectrum, leaving the emission line spectrum. The broad line region is then fit with a 
mixture of functions (such as multiple Gaussians or Gauss-Hermite polynomials). The continuum 
luminosity and line width are then extracted from the best-fit model. The measurement errors 
from the multiple component fits a re often estimated using some Monte Carlo methods (e.g., 
Shen et alJEoiU IShen & Liull2012h : mock spectra are generated either by adding noise to the 
original spectrum, or by adding "scrambled" residuals from the data minus best-fit model back to 
the model. The mock spectra are then fit with the same fitting procedure, and the formal errors 
are estimated from the distributions of the measured quantity from the mocks. This mock-based 
error estimation approach takes into account both the noise of the spectrum and ambiguities in 
decomposing different components in the fits. 



Empirical iron emission templates in the rest-frame UV to optical can be found in, e.g., 'Boroson & Greed <1992l) ; 
IVestergaard & W ilkes (2001); Tsuzuki et al. (2006). Using different iron templates may lead to small systematic offsets 
(< 0.05 dex) in the measured continuum luminosity and line width (e.g., Nobutaetal. 2012). Occasionally a Balmer 
continuum component is added in the pseudo-continuum fit to improve the fit around the "small blue bump" region near 
3000A (e.g., Grandi 1982; Dietiich et al. 2002), but s uch a component is generally difficult to constrain from spectra with 
limited wavelength coverage (see discussions in, e.g., IWang et alj|2009l ; IShen & Liil2012f) . 
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Below are some additional notes regarding continuum and line measurements. 



Narrow line subtraction Since the narrow line region (NLR) dynamics is not dominated 
by the central BH gravity, we want to subtract strong narrow line component before we 
measure the broad line width from the spectrum. This is particularly important for FWHM 
measurements, while for crune the effects of narrow lines are less important. For H/? and Ha, 
reliable constraints on the velocity and width of the narrow components can be obtained 
from the adjacent narrow hnes such as [O iii] /l/i4959,5007 and [S ii] /l/i67 17,6731. For 
Mg II and Civ, this is not so simple mainly for two reasons: 1) there are usually no adjacent 
strong narrow lines such as [O m] to provide constraints on the narrow line component; and 
even if [O iii] can be covered in other wavelengths there is no guarantee the NLR properties 
are the same for [O iii] and for Mg ii/Civ. 2) Although some quasars do show evidence of 
narrow component Mg n a nd Civ, it is unclear if this applies to the general quasar popu- 



lation. iShen & Liul (120121) found that for the 60 high-luminosity (L5100 > lO"*^ "* erg s~^) 
quasars in their sample with optical and near-IR spectroscopy covering Civ to [Om], the 
contribution of the narrow line component to Civ is too small to affect the estimated broad 
Civ FWHM significantly. However, for less luminous objects, the relative import ance of 



the na rrow line component to Civ might be larger (e.g., .Bachev et al...2004: .Sulentic et al 



20071). 



Remedy for absorption Sometimes there are absorption features superposed on the spec- 
trum, which is most relevant for Civ, and then Mg 11. Not accounting for these absorption 
features will bias the continuum and line measurements. While for narrow or moderately- 
broad abs orption trou ghs, manual or automatic treatments can greatly minimize their ef- 



fects (e.g., Shen et al.ll201 lb . there is no easy way to fit objects that are heavily absorbed 



(such as broad absorption line quasars). 

Effects of low signal-to-noise ratio (S/N) The quality of the continuum luminosity and line 
width measurements decreases as the quality of the spectrum degrades. In addition to 
increased measurement errors, l ow S/N data may also lead to biases in the spectral meas- 
urements. IPennev et al. fe009b) performed a detailed investigation on the effects of S/N on 



the measured HyS line width using many single-epoch spectra of two RM AGNs (NGC 5548 
and PG 1229 -1-204). They tested both direct measurements and Gauss-Hermite polynomial 
fits to the spectrum, and found that the best-fit line width is systematically underestimated 
at low S/N for both direct measurements and functional fits. The only exception is that 
their Gauss-Hermite fits to degraded NGC 5548 spectra tend to overestimate the FWHM 
at lower S/N. However, this is mainly caused by the fact that the Gauss-Hermite model is 
often unable to accurately fit the complex H/3 line profile of NGC 554 8. Using rnultiple 



Gaussian model fits, and for a much larger sample of SDSS quasars, IShen et al.l (1201 Ih 
also investigated the effects of S/N on the model fits by artificially degrading high S/N 
spectra (see their Figs. 5-8). They found that the exact magnitude of the bias depends on 
the line profile as well as the strength of the line. The continuum is usually unbiased as 
S/N decreases. The FWHMs and equivalent widths (EWs) are biased by less than ±20% 
for high-EW objects as S/N is reduced to as low as ~ 3/pixel. For low-EW objects, the 
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FWHMs and EWs are biased low/high by > 20% for S/N< 5/pixel. But the direction of 
the bias in FWHM is not always underestimation. 



3.3 Consequences of the Uncertainties in SE Mass Estimates 

Given the many physical and practical concerns discussed in § ^3. H and 13. 21 one immediately 
realizes that these mass estimates, especially those SE mass estimates, should be interpreted with 
great caution. Almost everyone acknowledges the large uncertainties associated with these mass 
estimates, but only very few are taking these uncertainties seriously. Since at present there is 
no way to know whether or not the extrapolation of these SE methods to high-z and/or high- 
luminosity quasars introduces significant biases, let us assume naively that these SE estimators 
provide unbiased mass estimates in the average sense, and focus on the statistical uncertainties 
(scatter) of these estimators. 

In mathematical terms, we have: 



me\m — m + G{0, ose) , (7) 

where m^, = log Mbh.se is the SE mass estimate, m = Mbh is the true BH mass, and G(ju, cr) is a 
Gaussian random deviate with mean ji and dispersion cr. I use x\y to denote a random value of x 
at fixed y drawn from the conditional probability distribution p{x\y). Eqn. ^ thus means that the 
distribution of SE masses given true BH mass, po{me\m), is a lognormal with mean equal to m 
and dispersion ose- It is then clear that this equation stipulates our assumption that the SE mass 
is on average an unbiased estimate of the true mass, but with a statistical scatter of ctse ~ 0.5 
(dex), i.e., the formal uncertainty of SE masses. 



3.3.1 The Malmquist-type bias (Eddington bias) 



Now let us assume that we have a mass-selected sample of objects with known true BH masses, 
and "observed" masses based on the SE estimators. By "mass-selected" I mean there is no se- 
lection bias caused by a flux (or luminosity) threshold - all BHs are observed regardless of their 
luminosity. If we further assume that the distribution of true BH masses in this sample is bottom- 
heavy, then a stati s tical bias in the SE masses naturally ari s es from the errors of SE masses (e.g., 
Shen et alJlIoOSal: iKellv eTal]|2009al: [shen & Kellvll2010l : iKellv et alj|2010l) . because there are 
more intrinsically lower-mass objects scattering into a SE mass bin due to errors than do in- 
trinsically higher-mass objects. This statistical bias can be shown analytically assuming simple 
analytical forms of the distribution of true BH masses. Suppose the underlying true mass distri- 
bution is a power-law, dN/dMBH ^^^^ Bayes's theorem tells us the distribution of true 
BH masses at given SE mass is (recall pQ{me\m) is the conditional probability distribution of nie 
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given to): 



I 



[m - (nie + ln(10)rMO"sg)]^ ' 



Thus the expectation value of true mass at given SE mass is: 

= TOe + ln(10)yMO"sE 



. (8) 



(9) 



Therefore for bottom-heavy (y^ < 0) true mass distributions, the average true mass at given 
SE mass is smaller by - ln(10)yMCsE '^^^ '■^^'^ '■^^ mass. This has an important consequence 
that the quasar black hole mass function (BHM F) constructed using SE virial masses will b e 
severely overestimated at the high-mass end (e.g. jKellv et al ] |2OO9all2O10HShen & Kellvll2012h . 



This statistical bias due to the uncertainty in the mass estimate s and a non-flat true mass distri 
bution is formally known as the Eddington bias ( Eddingtonll913 ). Historicall y this has also been 



referr ed to as the Malmquist bias in studies involving distance estimates (e.g., Lvnden-Bell et al 



1988h. which bear some resemblance to the familiar Malmquist bias in magnitude-limited samples 
(e.g., Malmquist 1922h. For th is re ason, this bias was ca lled the "Malmquist" or "Malmquist- 
type" bias in IShen et all dlOOSa") and lShen & KeUvl (1201(f). and I adopted this name here as well 
Perhaps a better name for this class of biases is the "Bayes correction", which then also ap- 
plies to the generalization of statistical biases caused by threshold data and correlation scatter (or 
measuremen t errors). The luminosity-dependent bias discussed next, and the Lauer et al. bias 
dLauer et al.ll2007i) discussed in ^4.31 can also be described by this name. 



3.3.2 Luminosity-dependent bias in SE virial BH masses 



Now let us take one step further, and consider the conditional probability distribution of to^ at 
fixed true mass to and fixed luminosity / = logL, p(me\m, 1). If the SE mass distribution at given 
true mass is independent on luminosity, then we have p{mi,\m, I) = p(me\m). This means that the 
SE mass is always unbiased in the mean regardless of luminosity. However, one may consider 
such a situation where pimdm, I) + p(me\m), which means the distribution of SE masses will be 
modified once one limits on luminosity. This is an important issue, since essentially all statistical 
quasar samples are flux-limited samples (except for heterogeneous samples, such as the local RM 
AGN sample), and frequently the SE mass distribution in finite luminosity bins is measured and 
interpreted. 

Below I will explore this possibility and its consequences in detail. To help the reader under- 
stand these issues, here is an outline of the discussion that follows: 1) I will first formulate the 
basic equations to understand the (mathematical) origin of the uncertainty in SE mass, ctse; 2) I 
will then provide physical considerations to justify this formulation; 3) The conditional probab- 
ility distribution of SE mass at fixed true mass and luminosity p(me\m, I) is then derived, and I 
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demonstrate the two most important consequences: the luminosity-dependent bias, and the nar- 
rower distribution of SE masses at fixed true mass and luminosity than the SE mass uncertainty 
o"se; 4) 1 then discuss current observational constraints on the luminosity-dependent bias and 
demonstrate its effect using a simulated flux-limited quasar sample. 

1} Understanding the origin of the uncertainty o"se in SE masses 

1 will use Gaussians (lognormal) to describe most distributions and neglect higher-order mo- 
ments, mainly because the current precision and our understanding of SE masses are not sufficient 
for more sophisticated modeling. Assuming the distributions of luminosity and line width at given 
true mass m both follow lognormal distributions, we can write such distributions as 

l\m = </>„ + Gi(0, a,), w\m = {w),„ + GiiQ, cr^) , (10) 

where notations are the same as in Eqn. (|7}, w = log W, and (),„ indicates the expectation value at 
m. The dispersions in luminosity and line width at this fixed true mass should be understood as 
due to both variations in single objects (i.e., variability) and object-by-object variance. The SE 
mass estimated using / and w are then (e.g., Eqn.|4]i: 

nielm - bl-\- cw -\- constant , (1 1) 

where the last term "constant" absorbs coefficient a and other constants from SE mass calibra- 
tions. Now let us consider the following two scenarios: 



A. In the ideal case where the SE method gives a perfect mass estimate, we have 

Gi(0,cr,) = -^G2(0,crJ, (12) 
b 

i.e., the deviations in luminosity and line width from their mean values at given true mass 
are perfectly correlated. The resulting distribution thus peaks at m with zero width, i.e., 
me\m — m — b{l}„, + c{w},„ -H constant. 

B. In the realistic case, certain amount of the deviations in / and w from their mean values are 
uncorrected with each other (either intrinsic or due to measurement errors). Without loss 
of generality, we can rewrite Eqn. (ITHt as 

l\m = <Z),„ + Gi(0, cr',) + Go(0, ct,,„), w\m = {w},„ + CM cr'J - ^Go(0, CT^.^r) , (13) 
where the total dispersions in the distributions of I and w are 



O"; = ycrj2 + crl^^ , cr^, = -^o-'^ + (o-^on-b / c)"^ ■ (14) 

Eqn. (fT3T l stipulates that some portions of the dispersions in I and w, described by cr^.a„, are 
correlated and do not contribute to the dispersion (scatter) of nig at m. On the other hand. 
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the remaining dispersions in / (cr^ and in w (cr^) are stochastic terms, and they combine to 
cause the dispersion of m^, at m: 

me\m = m + bGi(0, cr',) + cGiiO, cr'J^m + G(0, o-^e) , (15) 

where 

f^SE = ^{b(r',Y + (ccT'S (16) 
is the formal uncertainty in SE mass, i.e., the scatter in nie at given true mass m. 

Eqn. ( fTSl ) through ( fTSI l provide a general description of SE mass error budget from luminosity 
and line width, and form the basis of the following discussion. From now on I will only consider 
the realistic case B. 



2) Physical considerations on the variances cr'^, cr'^ and 0-^0^ 

Most of the studies to date h ave imphcitly assum ed cr ' = in Eqn. ([13]l. wi th the few excep- 
tions in e.g.. lShen et"al](l2008al) . IShen & KeilvlfeOlOl) and lShen & Kellvl (I2OI2I) . a-', = imposes 
a strong requirement that all the variations in luminosity are compensated by line width such that 
the uncertainty in SE masses now completely comes from the cr[^, part in line width dispersion. 
While this is what we hope for the S E method, there are ph ysical and practical reasons to expect 
a non-zero crj, as discussed in, e.g.. IShen & Kelly (l2012h . Specifically we have the following 
considerations: 



(a) the stochastic continuum luminosity variation and response of the BLR (hence the re- 
sponse in line width) are not synchronized, as resulting from the time lag in the reverberation of 
the BLR. The rms continuum variability on timescales of the BLR light-cross time is ~ 0.05 dex 
using the ensemble structure function in, e.g., MacLeod et al. (2010); 

(b) even with the same true mass, individual quasars have different BLR properties, and 
presumably the measured optical-UV continuum luminosity is not as tightly connected to the 
BLR as the ionizing luminosity. Both will lead to stochastic deviations of luminosity and line 
width from the perfect correlation (source-by-source variation in the virial coefficient /, scatter 
in the R - L relation, etc.). The level of this luminosity stochasticity is unknown but is at least 
0.2 - 0.3 dex given the scatter in the R- L relation alone, and thus it is a major contributor to crj; 

(c) although not explicitly specified in Eqn. ( fT3l l. there are uncorrelated measurement errors 
in luminosity and line width; typical measurement error in luminosity for SDSS spectra (with 
S/N~ 5 - 10/pixel, e.g., see fig. 4 of Shen et al. 2011) is ~ 0.02 dex (statistical only), but 
increases rapidly at low S/N; 

(d) and finally, what we measure as line width does not perfectly trace the virial velocity. This 
is a concern for essentially all three lines, and for both of the two common definitions of line width 
(FWHM and crune). Two particular concerns arise. First, single-epoch spectra do not provide a 
line width that describes the reverberating part of the line only, thus some portion of the line width 
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may not respond to luminosity variations. Second, if a line is affected by a non-virial component 
(say. Civ for instance), and if this component strengthens and widens when luminosity increases, 
the total line width would not response to the luminosity variation as expected. As in (b), this 
contribution to the uncompensated (by line width) luminosity variance crj is unknown, but could 
be as significant as in (b). 



One extreme of (d) would be that line width has nothing to do with the v irial velocity e xcept 
for providing a mean value in the calibrations of Eqn. (|4), as suggested bv ICroomI (1201 11) . i.e., 
CcoiT - 0. In this case while the average SE masses are unbiased by calibration, the luminosity- 
dependent SB mass bias at given true mass is maxi mum (see below). No te that this cr^on = case 
was already considered in I Shen et all (l2008al) and IShen & Kellvl (|2010|) when demonstrating the 
luminosity-dependent bias, and is only a special case of the above generalized formalism. 



On the other hand, a-^o„ > would mean that line width does respond to luminosity variations 
to some extent, justifying the inclusion of line width in Eqn. (|4]l. This was indeed seen at least for 
some local, low-luminosity objects, although not so much for the high-lumino sity SPSS sample, 
based on th e tests descr i bed in ^3.1.11 additional evidence is provided in, e.g. jKellv & Bechtold 
(l2007h and lAssef et al.l (120121) . again for the low-luminosity RM AGN sample. Therefore, the 
most realistic scenario is that at fixed true mass, some portions of the dispersions in luminosity 
(or equivalently, Eddington ratio) and in line width are correlated with each other, and they cancel 
out in the calculation of SE masses; the remaining portions of the dispersions in / and w are 
stochastic in nature and they combine to contribute to the SE mass uncertainty (as in Eqn.fTSIl. In 
other words, we expect crj > 0, cr^ > 0, and tTcorr > 0. For simplicity I take constant values for 
these scatters in the following discussion, but it is possible that they depend on true BH mass. 



3) The distribution of SE mass at fixed true mass and luminosity p(mg\m, l) 

Now that we have formulated the distributions of /, w and m^, at fixed m (e.g., Eqns. fTTIand 
[T3T l. we can derive the conditional probability distribution of nie at fixed m and I, p{mj\m, I). It is 
straightforward to showj^ (again using Hayes's theorem) that this distribution is also a Gaussian 
distribution, with mean and dispersion: 



{me}m,l ^m + 



cr, + O-t. 



-(/ - </),„) , 



{ba',fcTl 



• (co-'J 



(17) 



'^Here I give one possible derivation. For brevity I will drop m in all probability distributions, but it should be un- 
derstood that all these distributions are at fixed m. Consider the crj^ = case first, where we want to derive p(mg\l). 
Using Bayes's theorem, pinigll) oc p{mf)p(l\mg). We have pinig) oc /[2(*o-,) ] ({^^ all variance in nif comes 

from (tJ since cr'^ = 0), and p{l\mf) k e"''"' i> /(2o"cim) (j / ^-^^ only vary due to cTcon- at fixed ;«<,). Therefore 

p(me\l) oc e-<"''-<"«>"../)^/(2'^S), where (m,)„u is the same as in Eqn. (TT), and o-;^ = (hcr\fcrla,J{crf + o-2^„). Now add 
back in the <t'^, term, which will convolve p{mi\l) with a Gaussian distribution. Then the general distribution p{mi\l) for 
arbitrary values of it„, is also a Gaussian, with the same mean, but a dispersion that is broadened by ccrj^, (i.e., the same 
as in Eqn. ll7t . 
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Therefore we can generate the distribution of m^. at fixed m and / as: 

me\m, l^m+ /3(l - </>,„) + G(0, cr,„,) , (18) 
where (using Eqns.fl4land[T6Tl 

P = ,2^^' 2 = ^ ' = '^SE - = (b^'f + - P'oj . (19) 

O", + 0"corr O-/ 

Physically cr,„i is the dispersion of SE mass at fixed true mass and fixed luminosity. The parameter 
jS denotes the magnitude (slope) of the luminosity-dependent bias, and we hav^O < fi <b, where 
the lower and upper boundaries correspond to the two extreme cases crj = and (Xcon = 0. A 
larger jS means a stronger luminosity-dependent bias. Given the values of CTp cr'^ and (Tcon, and a 
SE calibration (b and c), all other quantities can be derived using Eqns. (fl4l i-(fT9t. 

Fig.|5]shows a demonstration with crj = 0.6, a-^o„ -Q.\,a-'^ -Q.\5,b - 0.5 and c = 2. In this 
example we have yS - 0.49, ctse = 0.42, and cr,„; = 0.3. The left two panels show the distributions 
of luminosity and line width at fixed true mass, from the stochastic term (crj, cr'^), the correlated 
term {cr^o„) and the total dispersion (cr/, cr^,). The right panel shows the distributions of SE masses 
at this fixed true mass for all luminosities (black dotted line) and for fixed luminosities (green and 
red dotted lines). The distributions of SE masses at fixed luminosity are both narrower and biased 
compared with the distribution without luminosity constraint. 

There are two important conclusions that can be drawn from Eqns. (fT7l l-(fT9li: 



• At fixed true mass m, the SE mass is over-/underestimated when luminosity / is higher/lower 
than the average value at fixed m. This is the "luminosity-dependent bias" first introduced in 
I shen et al. (2008a) and subsequently developed in Shen & Kelly (20i3) and[Shen & Kellvl 
(I2OI2I) . The magnitude of this bias, determined by fi, depends on how much of the dis- 
persion in luminosity at fixed true mass is stochastic with respect to line width (i.e., not 
compensated by responses in w). 

• Secondly, the variance in SE masses at fixed true mass and fixed luminosity, cr^^ is gener- 
ally smaller than the uncertainty of SE masses, cr^^, by an amount of {ficrif-. A simple way 
to understand this is that the uncertainty (variance) of SE mass at fixed true mass comes 
from the stochastic variance terms in both luminosity and line width. Therefore when one 
reduces the variance in either luminosity or line width by fixing either variable, the variance 
in SE mass is also reduced. This has important consequences in interpreting the observed 
distribution of SE masses for quasars in flux-limited samples or in narrow luminosity bins . 



aistribution or it, masses tor quasars m t lux-limited samples or m narrow luminosity bms . 
One cannot simply argue (as did in, e.g., iKollmeier et ai]l2006l:lsteinhardt & Elvisll2010bh 



that the uncertainty in SE masses is small because the distribution of SE masses for samples 



'For the sake of completeness, 1 note that fi > h could happen, if the line width were actually positively correlated to 
luminosity in the (T^on terms in Eqn. | |13) . Of course such a scenaiio is counter-intuitive (regarding the virial assumption) 
and thus unlilcely, and it means one should not use line width at all in estimating SE masses. 




Figure 5. Simulated distributions of / (left), w (middle), and (right) at fixed true mass m, following the 
description in 93.3.2l (e.g.. see Eqns. [T3ll9t . The example shown here assumes b = 0.5, c = 2 in Eqn. l fT3t . 
e.g., the typical values for SE mass estimators. Left: The distribution of luminosity / = logL. The black 
dotted line indicates the dispersion in crj = 0.6, the cyan dashed line indicates the dispersion in <Tcon = 0.1, 

and the black solid line indicates the total dispersion in o"/ = -^crj- + cr^^^ (which essentially overlaps 
with the dotted line given that ct; is dominated by crj). Middle: The distribution of line width w = log W. 
The black dotted line indicates the dispersion in cr'^, = 0.15, the cyan dashed line indicates the dispersion in 
-Cconb/c (i.e., the part that correlates with luminosity), and the black solid line indicates the total dispersion 
in CTiu = ^jcr'^ + (cr^anb/c)-. Again the solid line and the dotted line are almost on top of each other. Right: 
The distribution of SE mass m^, = log Mse. The black dotted line indicates the total dispersion of at fixed 
m, (TsE = .^(bo-'iy + (c(T^)- = 0.42. The green (red) dotted line indicates the distribution of at fixed 
m and fixed / = - 0.5 ((/)„, + 0.5), which is a Gaussian described by Eqns. ( I18t -J19t. In this example 
most of the dispersions in luminosity and line width are uncoiTelated with each other, leading to a large 
uncertainty in SE masses ctse = 0.42 dex. The inferred luminosity-dependent bias has a slope of /3 = 0.49. 
The dispersion of at fixed luminosity is only cr,„/ = 0.3 dex (Eqn.|19b, much smaller than the uncertainty 
in SE masses, ctse- 



with restricted luminosity ranges is narrow. The example shown in Fig. |5]clearly demon- 
strates that one can easily get a much narrower SE mass distribution at fixed true mass and 
luminosity, than the nominal SE mass uncertainty ose- On the other hand, if enforcing 
cr'i - 0, then ose = co"^ < ccr„„ and there will be tension betwe en the observed narrow dis- 
tribution in line width (e.g.,|Shenetal ]|2008al: l iFine et al ]|2008h . which indicates cr„, < 0.15 
dex for Mg ii, and the expectation that ose > 0.3 dex. 



4) Observational constraints on the luminosity-dependent bias 



The exact value of f} is difficult to determine observationally, although some rough estimates 
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Figure 6. A test on a non-zero fH using RM data for NGC 5548. Plotted here is the dependence of the 
virial products computed from 5100A continuum luminosity and line width as a function of luminosity, for 
NGC 5548 and for Hj8 only, for FWHM (left) and (Tii„e (right), respectively. This form of the virial product 
represents SE estimators with b = 0.5 and c = 2 in Eqn. and I am using luminosity instead of time lag t 
in computing the virial product since I am testing SE mass estimators. The measurements were taken from 
ICollin et al 1 l l2006l) . which are based on both mean and rms spectra during each monitoring period. Error 
bars represent measurement errors. The error bars in luminosity have been omitted in the plot for clarity. 
The continuum luminosity has been corrected for host starlight using the correction provided bv lBentz et al.l 
The black and blue dashed lines are the best linear-regression fits using the Bayesian method 
of lKellvl llOOi). for measurements based on mean and rms spectra, respectively. The residual correlation 
between the virial product and luminosity cannot be completely removed, and a positive p ~ 0.2 - 0.6 is 
inferred in all cases, although the uncertainty in j8 is too large to rule out a zero /? at > Scr significance. 
Figure adapted from lShen & Kellvl ( |2012|) . 



can be made based on monitoring data of single objects, or samples of AGNs with known "tiTie" 
mass (using RM masses or Mbh - o", masses). The former test constrains the stochasticity in 
single objects, while the latter test explores object- by-object variance. U sing the intensively 
monitored RM data for a single object, NGC 5548. IShen & Kellvl ( l2012h tested the possibility 
of a non-zero /3. The continuum luminosity of this object varied by ~ 0.5 dex within a decade, 
providing a test on how good line width varies in accordance to luminosity variations for a single 
object and for a single line. Fig. |6] shows the change of SE masses as a function of the mean 
continuum luminosity in each monitoring year, computed using both FWHM and o-\i„s from both 
mean and rms spectra in each year. There is an average trend of increasing the SE masses as 
luminosity increases in all four cases, although the trend is le ss obvious fo r criine-based SE masses. 
The inferred value of p, using the linear regression method in lKellvl (l2007'). is ~ 0.2—0.6, although 
the uncertainty in /3 is generally too large to rule out a zero yS at > 3cr significance. 
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Figure 7. Tests of a non-zero /3 using two-epoch spectroscopy from SDSS (She, Shen, et al., in preparation). 
The data and notations are the same as in Fig.[T] The SE masses remain more or less constant as luminosity 
changes only for the low-luminosity and low-z sample based on H/3. For the high-luminosity and high-z 
samples based on Mg ii and CiV, a luminosity-dependent bias in SE masses is inferred, with an error slope 
ofj8~ 0.5. 
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A similar test is based on the repeated spectroscopy in SDSS (see discussion in ^3.1. 11 1. 
While most objects do not have a large dynamic range in luminosity variations in two epochs, the 
large number of objects allows a reasonable determination of the average trend of SE masses with 
luminosity, for the whole population of quasars. In addition, we want to include measurement 
uncertainties (both statistical and systematic), which allows us to make realistic constraints, as 
measurement errors will always be present. As shown in Fig. [T] the line width from single-epoch 
spectra does not seem to respond to luminosity variations except for low-luminosity objects based 
on Hy6, as expected from the physical/practical reasons I described above. I plot the changes in SE 
masses as a function of luminosity changes in Fig. Q From this figure I estimate p ~ 0.5 for the 
high-L samples based on Mg ii and Civ, and j3 ~ for the low-L sample based on H/3 (She, Shen 
et al., in preparation). The difference in the low-L (low-z) and high-L (high-z) samples could be 
due to a luminosity effect, e.g., the correspondence between line width and luminosity variations 
is poorer at higher luminosities, or due to the difference between H/3 and the other two lines (She, 
Shen et al., in preparation). 



Shen & Ke^^^ ll2012t) also attempted to constrain yS using forward Bayesian modeling of 
SDSS quasars in the mass-luminosity plane (see ^4.2b . While the results suggested a non-zero 
/? ~ 0.2 - 0.4, the constraints were not very strong (see their fig. 1 1). Combining all these tests, 
we can conclude the following: /3 is probably smaller than ~ 0.5 (i.e., line width still plays some 
physical role in SE mass estimates) but unlikely zero, although the exact value is uncertain. The 
value of yS likely also depends on the specific line. More monitoring data of individual AGNs, 
and/or a substantially larger sample of AGNs with RM mass (or Mbh - cr« masses) spreading 
enough dynamic range in luminosity at fixed mass, will be critical in better constraining /3. 

The effec ts of the lu minosity-dependent bias on flux-limited samples are dis cussed intens- 



ively in, e.g., IShen et all QOOSa), Shen & Kelly (2010) and Shen & Kell^ (l2012h . Here I use a 



simple simulation of a power-law true BH mass distribution, dN/dMsH ^mf' demonstrate 
these effects in Fig. |8] This steep true mass distribution was chosen to reproduce the distribution 
at the high-mass end for SDSS quasars (S hen et al. 2008a), which is certainly not appropriate at 
the low-mass end. I use the same example as in Fig. |5] for all the dispersion terms in luminosity 
and line width at fixed true mass. The true masses are distributed between 5x10'' M© and 10'" Mq 
according to the specified power-law distribution. The mean luminosity at fixed m is determined 
assuming a constant Eddington ratio A = 0.05. Then the instantaneous luminosity and SE mass 
at each true mass m are generated using Eqns. (fTJt-(fT9). The resulting distribution in the m- nie 
plan is shown in black contours and points in Fig. |8] where I also show the distribution of a flux- 
limited (luminosity-limited) subset of BHs with / > (/),„=8.3, the mean luminosity corresponding 
to Mbh = 2 X 10^ Mq. The simulated distributions in luminosity, line width and SE virial masses 
are consistent with those for SDSS quasars when a similar flux limit is imposed on the simulated 
BHs. It is clear from Fig.[8]that for the flux-limited subset, the SE masses are biased high from 
the true masses. This is because in this simulation we have (3 - 0.49, which implies a significant 
luminosity-dependent bias at fixed true mass. Then the bottom-heavy BH mass distribution and 
the large scatter of luminosity at fixed true mass lead to more overestimated SE masses scattered 
upward than underestimated ones scattered downward, causing a net sample bias in SE masses 



up 

Isi 



Shen et al.,.2008a:.Shen & Kellv..20ia,2012) . 
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Figure 8. A simulated population of BHs with true masses Mbh within 5 x 10' - 10'" M©, distributed as a 
power-law dN/dMBH °^ ^Bif- I have assumed a constant Eddington ratio A = 0.05 to generate the mean 
luminosity (/),„ at fixed true mass. I then used the dispersions specified by the example shown in Fig.|5]to 
generate instantaneous luminosity and SE masses. The resulting distributions in luminosity, line width and 
SE masses are consistent with those for the SDSS quasar sample (e.g., IShen et alj2oTlh . The black contours 
(local point density contours) and points show the distribution in the - m plane. The green line shows 
the unity relation lUg = in. As expected, there is substantial scatter in irig at fixed m due to the uncertainty in 
SE masses (ctse = 0.42 dex). The red contours and points show the distribution for a subset of quasars with 
/ > (/>„,=8.3, the corresponding mean luminosity at Mbh = 2 x 10** Mq (marked by the green circle). The SE 
masses are biased high from their true masses for this flux-limited sample (i.e., most points are above the 
unity relation) due to the substantial luminosity-dependent bias, the large dispersion in luminosity at fixed 
true mass, and the bottom-heavy true mass distribution in this example. 
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4. Applications to Statistical Quasar Samples 

Despite the many caveats of SE mass estimators discussed above, they have been extensively 
used in recent years to measure BH masses in quasar and AGN samples over wide luminosity 
and redshift ranges, given their easiness to use. These applications include the Eddington ratio 
distributions of quasars, the demographics of quasars in terms of the black hole mass function 
(BHMF), the correlations between BH mass and host properties, and BH mass dependence of 
quasar properties. It is important to recognize, however, that these BH mass estimates are not true 
masses, and the uncertainty in these mass estimates has dramatic influences on the interpretation 
of these measurements. 

Below I discuss several major applications of the SE virial mass estimators to statistical 
quasar samples. Other applications of these SE masses, such as quasar phenomenology, while 
equally important, will not be covered here. 



4.1 Early Growth of SMBHs 



One of the strong drivers for developing the SE virial mass technique is to estimate BH masses 
for high redshift quasars to better than a factor of ten accuracy, and to study the growth of SMBHs 
up to very high redshift (e. s.. IVestergaardl2004l) . Such investigations have been greatly improved 
in the era of modem, large-scale spectroscopic surveys. The SDSS survey has been influential 
on this topic by providing more than tens of thousands of optical quasar spectra and SE mass 
estimates up to z ~ 5 (e.g., McLure & Dunlopll2004i iNetzer & Trakhtenbrotll2007[ IShen et al 



2008all201 lULabita et alJl2009allbh . On the other hand, deeper and dedicated optical and near-IR 



spect r oscopic prograrns are probing the SMBH growth to even higher redshift (e.g., iJiang et al 
20071: iKurk et al.ll2007l: IWiUott et alJl2010t iMortlock et all201 ih . 



In Fig.|9]l show a compilation of SE virial mass estimates for quasars over a wide redshift 
range < z < 7 from diflrerent studies. The black dots show the SE masses from the SDSS DR7 
quasar sample dSchneider et al. 2010l; Shen et al. 201 1 ), which were estimated based on H/? for 
z < 0.7, Mgii for 0.7 < z < 1.9 and Civ for 1.9 < z < 5. As discussed in ^3.1.71 the rehabihty of 
Civ-based SE masses for the high-z and high-luminosity quasars has been questioned, and several 
studies have obtained near-IR spectra for z > 2 quasars to get HyS-based (filled symbols) and 
Mg ii-based (open symbols) SE masses at high redshift. Albeit with considerable uncertainties 
and possible biases in individual SE mass estimates (all on the order of a factor of a few), these 
studies show that massive, > 10^ Mq BHs are probably already in place by z ~ 7, when the age of 
the Universe is only less than 1 Gyr The abundance of these massive and active BHs then evolves 
strongly with redshift, showing the rise and fall of the bright quasar population with cosmic time. 



One outstanding question regarding the observed earliest quasars is how they could have 
grown such massive SMBHs given the limite d time they have, which is a non-triv ial problem 
since the discoveries of z > 4 quasars (e.g.. lTurnerlll99U iHaiman & Loebll2001l) . One con- 
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cem is if these highest redshift quasars have their luminosities magnified by gravitational lensing 
or if their luminosities are strongly beamed (e.g., Wyithe & Loeb 2002; Haiman & Cen 200^), 
which will affect earlier estimates of their BH masses using the Eddington-limit argument. The 
lensing hypothesis will also lead to overestimated virial BH masses. However later deep, high- 
resolution imaging of z > 4 quasars with HST did not find any multiple images around these ob- 
jects (e.g..iRichards et al.ll2004bi l2006bl) . rendering the lensing hypothesis highly unlikely (e.g., 
Keeton et alJl200 5S. Strong beaming can also be ruled out bas ed on the h igh values of the ob- 



served line/continuum ratio of these high-redshift quasars (e.g.. lHaiman & Cen 2002) 



Given the e-folding time introduced in ^ = 4.5 x 10** ^^^-j^ yr, and a seed BH mass Mseed 
at an earlier epoch z,, the final mass at z/ ~ 6 is 



A^BH = Mseed exp 



tf - u 



(20) 



where tf and f, are the cosmic age at z/ and z,, respectively. 



Assuming continuous accretion with constant radiative efficiency e and luminosity Eddington 
ratio A and without mergers, I showed in Fig. |9]three different growth histories from a seed BH 
at higher redshift. The solid lines are for a seed BH at z = 20 and A - \, i.e., Eddington- 
limited accretion; the dashed lines are for a seed BH at z = 30 and A - \\ the dotted lines are 
for a seed BH at z = 20 and A - 1.5, i.e., mildly super-Eddington accretion. For each model 
I used three seed BH mass, Mseed = 10,20, lOOM©, which encloses the reasonable ranges of 
predicted rernnant B H mass from the first generation of stars (Pop III stars, for a review see, e.g., 
Bromm et al.ll2009h . Then it is clear from Fig. |9]that, if the accretion is Eddington-limited, it is 
difficult to grow > 10^ BHs at z ~ 6 from a Pop III remnant seed BH at z ~ 20 - 30 (where 
such first stars were formed out of ~ 10^ Mq halos, corresponding to ~ 3 - 4cr peaks in the 
density perturbation field). On the other hand, if allowing mildly super-Eddington accretion, then 
a ~ 10"^ Mq BH can be readily formed at z ~ 7 from a large Pop III star remnant Mseed ~ 100 Mq, 
although more recent simulations suggest somewhat lower mas ses of Pop III stars due to possible 
effects of clump fragme ntation and/or radiative feedback (e.g.. lTurk et al ] l2009t iHosokawa et al 



2011l : IStacv et al . 2012), and hence a lower typical value of the remnant mass of less than tens 
of solar masses. Mildly super-Eddington accretion (up \a A ~ few) could happen, for instance, 
if the radiation and density fields of the accretion flow are anisotropic and most of the accretion 
flow is not impeded by the radiation force. Mergers between BHs at high-z can also help with the 
required growth if the coalesced BH is not ejected from the halo by the gravitational recoil from 
the merger The main challenge here is whe ther or not such critica l accretion can maintain stable 
and uninterrupted for the entire time (e.g., iPelupessv et al. 112007^. But in any case, it is quite 
likely that the observed z > 6 quasars are all bom in rare environments of the early Universe, 
thus extreme conditions (such as large gas density, high merger rate, etc.) may have facilitated 
their growth. Indeed, some theoretical studies can successfully produce such massive BHs at 
z ~ 6 growing from a typical Pop III remnant BH seed without super-Eddington accretion (e.g., 
Yoo & Miralda-Escudell2004t iLi et al.ll2007t iTanaka & Haiman| [2009). But the detailed physics 
(accretion rate, mergers, BH recoils, etc.) regarding the formation of these earUest SMBHs is still 
uncertain to some large extent. 
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Figure 9. A compilation of SE virial mass estimates from different samples of quasars. The black dots are 
the SE masses for SDSS quasars from IShen etaP (1201 ih . based on H/3 (z < 0.7), Mgii (0.7 < z < 1.9) 
and Civ (z > 1.9). Given the potential caveats of CiV-based SE masses, near-IR spectroscopy has been 
undertaken to estimate SE masses for z > 2 quasars based on Hj8 and Mg ii. The different large symbols 
are for SE masses based on H/3 (filled symbols) and Mg ii (open symbols) from lShemmer et all ( l2004i filled 
squares), Netzer et al. ('2007', filled circiles), 'Trakhtenbrot et al. ' ('20 ll", open circles), iKurk et al] ( l2007l 
open triangles), Willott et al, ( 201(1 open squares), and Mo rtlock et al.. (.201 L open circle with cross). Bear 
in mind the large uncertainty associated with individual SE masses and potential biases. I also show the 
predicted SMBH growth at z > 6 based on simple, continuous accretion models with constant Eddington 
ratio A and radiative efficiency 6 = 0.1, as described in Eqn. ( I20t . The solid lines are Eddington-limited 
(A = 1) growth models from a seed BH at z = 20; the dashed lines are A = I growth models from a seed 
BH at z = 30; the dotted lines are mildly super-Eddington (A = 1.5) growth models from a z = 20 seed BH. 
For each model I used three seed BH masses, Mseed = 10, 20, 100 M© to accommodate reasonable ranges of 
seed BH mass from a Pop III star remnant at z ~ 20 - 30. 
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While this is not seen as an immediate crisis, there are multiple pathways to make it much 
easier to grow > 10^ Mq SMBHs at z > 6 by boosting either the accretion rate or the seed BH 
mass (for a recent review, see, e.g., IWonterili2010[ lHaimanll2012l) . These recipes include: 



1) supercritical accretion (e.g., Volonteri & Reesll2005 1 where the accretion rate Mbh greatly 
exceeds the Eddington limit with a canonical ra diative efficiency e - 0. 1 . One possibih ty is that 
the radiation is trapped in the accretion flow (e.g., lBege lmanll979l:IWvithe & Loebl2012h . leading 
to a very low e and hence a much shorter e-folding time. Note that in such a radiatively inefficient 
accretion flow (RIAF), the luminosity is still bounded by the Eddington limit; 

2) rapid formati on of massive (~ 10^ - 10^ Mp) BH se eds from direct colla pse of primor- 
dial gas clouds (e.g., Bromm & Loebl 2003 : Begelman et al. 2006; Aaarwal et al. 2012) or from 



a hypothetical supermassi ve star or "quasi-star" (e.g., IShibata & Shapiro 2002: Begelman et al 



2008 : Johnson et al.ll2012 ) at high redshift. Supercritical accretion may also be expected in some 
of these models to grow to the final seed BH mass, which then continue to accrete in the normal 
way. By increasing Mseed it requires much less e-folds to grow to a > 10'^ BH. Another pos- 
sible route to produce massive BH seeds up to ~ 10^ Mp is by the runaway coUisional growth in 
a dense star cluster formed in a high-redshift halo (e.g.. lOmukai et al ] l2008h . 



4.2 Quasar Demographics in the Mass-Luminosity Plane 



BH mass estimates provide an additional dimension in the physical properties of quasars. The dis- 
tribution of quasars in the two-dimensional BH mass-luminosity (M-L) plane conveys important 
information about the accretion process of these active SMBHs. The first quasar mass-luminosity 
plane plot was made by Dibai in the 1970s as mentioned in ^2.21 Over the years, such a 2D plot 
has been repeatedly generated based on increasingly larger quasar samples and improved BH 



nas been repeateaiy generateu paseu on increasmgly larger quasar samples ana irnproved an 
mass estimates (e.g.. Wandel et al.lll999t IWoo & Urrvlliool Kollmeier et al.ll2006l lShen et al 



2008aHShen & Kellvll2 012'). and the much better statistics now allows a more detailed and deeper 



look into this quasar mass-luminosity plane. 



In what follows I will mainly focus on the SDSS quasar sam ple because this is the largest and 
most homogeneous quasar samples to date. But as emphasized in lShen & Kellvl (120121) . the SDSS 
sample only probes the bright-end of the quasar population, and to probe the mass and accretion 
rate of the b ulk of qu asars it is neces sary to as semble deeper spectroscopic quasar s amples (e.g., 
KoUmeier et al. 2006: lGavignaud et al. 2008; Trump et al.lE009l:lNobuta et alJl2012l) . 



Since I have emphasized the distinction between true BH masses and SE mass estimates, I 
shall use the term "observed" or "measured" to refer to distributions based on SE mass estimates, 
to distinguish them from "true" distributions. Fig. [TOl shows such an observed mass-luminosity 
plane from the same collection of quasars as shown in Fig. |9] Note that these samples are flux- 
limited to different magnitudes, and several high-redshift samples based on H/3 or Mgii (i.e., 
large symbols) have a higher flux-limit than the SDSS. I also used slightly different values of 
bolometric corrections to convert continuum luminosity to bolometric luminosity for those non- 
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SDSS samples. From this plot we see that the observed distributions of quasars are bounded 
between constant Eddington ratios 0.01 < A < 1, with median values of {A} ~ 0.1 - 0.3 for 
SDSS quasars, and somewhat higher values for the z > 5 samples. The dispersion in Eddington 
ratio in these flux-limited samples is typically ~ 0.3 dex. Similar distributions were obse rved by. 



e.g., Kollmeieret al. (200^. However, as demonstrated in, e.g., Shen et al. (2008a); Kel lv et al. 
(l2009all2010h : IShen & Keilvl(l2012l) : [Kellv & Shed jlOlA the observed distribution suffers from 
the sample flux limit such that low-Eddginton ratio objects have a lower probability being selected 
into the sample, and from the uncertainties and statistical biases of SE masses relative to true 
masses. The selection effect due to the flux limit and errors in SE masses dramatically modify 
the intrinsic distribution in the mass-luminosity plane, and must be taken into account when 
interpreting the observations. 



The best approach to tackle these issues is a forward modeling, in which one specifies an 
underlying distribution of true masses and luminosities and map to the observed mass-luminosity 
plane by imposing the flux limit and relati ons between SE virial masses and true masses (e.g., 
Shen et al ]|2008al:lKellv et alJl2009aLEoiol) . Then the comparisons between model and observed 



distributions constrain the model parameters with standard Markov Chain Monte Carlo (MCMC) 
techniques and Bayesian inference. T his is a complicated and model-dependent pro blem, and 
the best efforts so far are the studies bvlShen & Ke 



2009a') and'K ellv et alj \20ld ). Alternatively, 
Schulze & Wisotz ki (2010i) developed a maximum likelihood method (also a forward modeling 



on earlier wo rk by Shen et alj (l2008ab . iKellv et al.l 



ivl(l2012l) a ndlKellv & Shen j2013h . building 



method), which accounts for the effect of the flux limit, but not the errors in SE masses (al- 
though the SE errors can be incorporated in such a framework as well). This maximum likelihood 
method was subsequently adopted in Npbuta et al. (2 01Z) when modeling a faint quasar sample 
(again, SE errors not taken into account). Most othe r quasar mass demographic studies, however . 



did not explic itly model either of these effects (e.g.. lGreene & 



iss demogi 
Hol l2007t 



Vestergaard et al.ll2008 



Vestergaard & Osmer.2009,) . 



Shen & Kellvl 1 120 12n used forward modeling with Bayesian inference to model the observed 



distribution in the mass-luminosity plane of SDSS quasars, taking into account a possible luminosity- 
dependent bias (i.e., P + Q, see ^3.3.2l i to be constrained by the data. The flux li mit of the SDSS 



sample is taken into account using publi shed se l ection functions of SDSS quasars (iRichards et al 



2006al) . Based on this approach. IShen & Kellvl (120 12h found e vidence for a non-zer o yS, although 



the constraints on p are weak and cannot rule out a null value. iKellv & ShenI (120131) used a more 
flexible model parametrization to describe the underlying true distr ibutions (in BH mass and Ed- 
dington ratio), but fixed /? = to test how sensitively the results in lShen & Kellvl (120 12h depend 
on different model parameterizations. They found that the main conclusions are generally con- 
sistent, although the results in the latter work are less constrained than in the former, due to the 
more flexible models. Both studies revealed that based on the SDSS quasar sample alone, it is 
difficult to constrain the BHMF to better than a factor of a few at most redshifts. This is both 
because the SDSS sample only probes the tip of the active SMBH population at high-z, and to a 
larger extent, because the errors of SE masses are poorly understood. However, there are some 
solid conclusions from the two studies: 
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Figure 10. The observed quasar mass-luminosity plane based on SE masses for quasars in a wide range of 
redshifts (0 < z < 7) from the samples compiled in Fig.|9] The dots are from the SDSS sample in.Shen et alj 
for quasars at z < 0.7 (H/3-based SE masses; green), 0.7 < z < 1.9 (Mgn-based SE masses; cyan), 
and 1.9 < z < 5 (CiV-based SE masses; red). The large symbols are from various z <: 2 samples using Hj8 
or Mg n-based SE masses. I have used slightly different bolometric corrections for these non-SDSS samples 
from those used in the original papers. The SDSS quasars have Eddington ratios 0.01 < A < I with a mean 
value of (A) ~ 0.1 - 0.3, while the other high-z samples have even higher Eddington ratios. Given the flux- 
limited nature of all these samples and the errors in SE masses, the observed Eddington ratio distribution is 
highly biased relative to the intrinsic distribution (see g4.2l for details). 
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Figure 11. The simulated mass-luminosity plane at z = 0.6 based on the modeling of SDSS quasars in 
|shen^feii3 iloil), which extends below the flux limit (the black horizontal line). The y-axis plots the 
restframe 2500A monochromatic luminosity, and the bolometric luminosity is Lboi ~ 5L. The red contour 
is the distribution based on true BH masses and is determined by the model BHMF and Eddington ratio 
distribution in IShen&Kellvl ( l2012l ). The black contour is the distribution based on H/3 SE BH masses. The 
flux limit only selects the most luminous objects into the SDSS sample (which are closer to the Eddington 
limit), and the distribution based on SE virial BH masses is flatter than the one based on true masses due to 
both the scatter cr,,,; and a non-zero /3 ~ 0.2 for this redshift bin (see ^Tlll and lShen & Kellvll2012h . 



1 . The observed distribution in the mass-luminosity plane is quite different from the intrinsic 
distribution, due to the flux limit and uncertainties in SE masses. This is demonstrated in 
Fig. [TTl w hich shows a quasar M - L plane at z = 0.6 based on the modeling of SDSS 
quasars by iShen & Kellvl (120121) . In this plot luminosity L is the restframe 2500A mono- 
chromatic luminosity, and the bolometric luminosity is Lboi ~ 5L. The red contours are the 
true distribution of quasars, while the black contours are the measured distribution based 
on HjS SE virial masses. The black horizontal line indicates the flux limit of the sample, 
hence only objects above this line would be selected in the SDSS sample, which form the 
observed distribution. The flux limit only selects the most luminous objects into the SDSS 
sample, missing the bulk of low Eddington ratio objects; even the highest mass bins are 
incomplete due to the flux limit. The distribution based on SE virial BH masses is flatter 
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than the one based on true masses due to the scatter and luminosity-dependent bias of these 
SE masses. 

For the observed distribution based on SE masses, there are fewer objects towards larger 
SE masses and luminosity. This was interpreted as the lack of massive black holes ac- 
creting at high Eddington ratios, or the so-called "sub-Eddington boundary" claimed by 
Steinhardt & ElvisI (l2010ah . However, such a feature is caused by the flux limit and errors 
in SE masses, and there is no evidence that high-mass quasars on average accrete at lower 
Eddington ratios, not for broad-line objects at leas0 This conclusion seems to be robust 
against different model p arameterizations of the underlying true distributions in the for- 
ward modeling approach (iKellv & ShenlEoTsI) . Of course here I am assumi ng no system- 
atic biases in these FWHM-based SE masses measured in iShen et al.l (1201 lb . It is possible 
that cTiine-based SE masses are more reliable, in which case there would be a "rotation" in 
the mass-luminosity plane using crjing-based SE masses, as discussed in ^3.1.3[ This also 
tends t o reduce this "sub-Eddington boundary" in the observed plane (e.g.. Rafiee & Hall 
2011allbh . 5ut a full modeling taking into account both the flux limit and SE mass errors is 
yet to be performed with crjine-based SE masses (i.e., the interpretation by Rafiee & Hall is 
still based on "observed" rather than true distributions). 



The intrinsic Eddington ratio distribution at fixed true mass is broader (~ 0.4 dex) than the 
oZ75erve(i Eddington ratio distribution in flux-limited samples (< 0.3 dex), and the mean Ed- 
dington ratio in the flux-limited samples based on SE masses is higheo than the mean Ed- 
dington ratio for a ll active SM BHs (most of which are not detected). This is consistent with 
earlier studies bv IShen et al.l (i2008a) and iKellv et al. (2010). Some de eper spectroscopic 



surveys indeed start to find these lower E ddington ratio objects (e.g.. iBabic et al.l 12007 



Gavignaud et al . "2008'; Nobuta et al. "2012'), and are consiste nt with the model-extrapolated 



distributions from Shen & Kellv (2012) and Kellv & Shen (120131) : however, since in gen- 
eral these deep data are noisier and the selection function is less well quantified than SDSS, 
care must be paid when inferring the dispersion in Eddington ratios for these faint quasars. 



The next step to utilize this quasar mass-luminosity plane is to measure the abundance of 
quasars in this plane, and study its redshift evolution. This is a much more powerful way to study 
the cosmic evolution of quasars than traditional ID distribution functions such as the luminosity 
function (LF) and the quasar BHMF 

I demonstrate the power of the mass-luminosity plane in quasar demographic studies in Fig. 
[T2I This^is the same simula ted, true quasar M - L plane at z = 0.6 as in Fig.[TT] using the mod- 
els in iShen & Kellvl (l2012h constrained using SDSS quasars. The 2D density (i.e., abundance) 



There could be a significant population of high-mass, non-broad-line AGNs accreting at low Eddington ratios, pos- 
sibly via a different accretion mode than broad line quasars. 

"These SE masses are on average overestimated due to the luminosity-dependent bias discussed in j|3.3.2l which 
tends to underestimated the true Eddington ratios. But the mean observed Eddington ratio based on SE masses of the 
flux-limited sample is still higher than the mean value for all quasars extending below the flux hmit (see fig. 19 of 
IShen & Kellvll2012h . 
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Fi gure 12. An examp le of the forward modeling of quasar demographics in the mass-luminosity plane 
bv lshen & Kdh] llOli) . modeled at z = 0.6. Left: The simulated mass-luminosity plane (with true BH 
masses), which extends below the SDSS flux limit (the black horizontal line). Shown here is the comoving 
number density map [0(Mbh, ^')], where only regions with<^(MnH,Q > 10""Mpc"MogL-'logMB;, 
are shown. The two diagonal lines indicate constant Eddington ratios of 1 and 0.01. The flux limit only 
selects the most luminous objects into the SDSS sample. Right: Projections of the 2D distribution onto 
the luminosity and BH mass axes, i.e., the LF (upper right) and BHMF (bottom right) of quasars at 
z = 0.6. The points are binned observational data and the lines are best-fit models: the magenta lines 
show the results for all broad-line quasars (corrected for the flux limit) and the green lines show those for 
the flux-limited sample. The thickness of the lines indicates the Icr model uncertainty. I have used the 
/-band absolute magnitude instead of bolometric luminosity in presenting the LF, and I have also taken into 
account the difference between true BH masses (colored lines) and SE virial BH masses (points). This 
forward-modeling framework accounts for the selection incompleteness in BH mass due to the flux limit, 
and the uncertainties in SE virial BH mass estimates, and can constrain the 2D distribution down to ~ 3 
magnitudes fainter than the flux limit dShen & Kellvll2012l) . However, to make more robust constraints at 
the faint end, deeper quasar samples are highly desirable. 

The 2D distribution of quasars offers significantly more information on the evolution of quasars than 
ID distributions (i.e., LF or BHMF) alone. Specifically if we focus on snapshots of the mass-luminosity 
plane in short time intervals (e.g., shorter than typical quasar lifetime), then objects will move along certain 
evolutionary tracks in the plane, determined by the evolution of Eddington ratio (possibly also coupled with 
the change of the radiative efficiency as accretion rate evolves). This evolution in the mass-luminosity plane 
can be modeled in great detail, and compared with either analytical models or numerical simulations of 
SMBH growth. This is somewhat similar to the simple arguments in the quasar LF literature (such as the 
"pure luminosity evolution" model), but the diagnostic capabilities of the mass-luminosity plane would be 
far more powerful and self-consistent in terms of SMBH growth and light curve evolution. 
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of quasars in this plane is shown in the color-coded contours. The traditional LF and BHMF 
(shown in the right panels) are then just the projection onto each axis. In the right panels I also 
demonstrate the differences between using true BH masses and SE virial masses, as well as the 
effect of the sample flux limit. It is clear from this demonstration that the ID distribution func- 
tions lose information by collapsing on one dimension, and a better way to study the demography 
of quasars is to measure their abundance in 2D, since the mass and luminosity of a quasar are 
physically connected by the Eddington ratio. The ultimate goal is to study the evolution of the 
quasar density in the M — L plane as a function of time. R ecent studies have started to work in this 
direction (e.g., IShen & Kellvll2012l;lKellv & Sheij 120131) . although deeper quasar samples and a 
better understanding of SE mass errors are needed to utilize the full power of the M - L plane. 



To summarize, the quasar mass-luminosity plane has great potential in st udying quasar evol- 



ution, an d efforts have been underway to investigate this plane in detail (e.g.,|St einhardt & Elvid 
2010allbl : lsternhardt et al.l2oTTllSteinhardt & Elvisl20l"TllShen & Kellvl2012l: lKenv & Shen 261J). 
However, it should always be kept in mind that the "observed" distribution in the M - L plane is 
not the true distribution. I strongly discourage direct interpretations of the observed distributions 
based on SE masses and flux-limited data, which can easily lead to superficial or even spurious 
results. 



4.3 Evolution of BH-Bulge Scaling Relations 

Another important application of SE virial mass estimators is to study the MfiH-host scaling rela- 
tions in broad line AGNs, and to probe the evolution of these relations at high redshift. Measuring 
the MflH-host relatio ns in low redshift quasars and AGNs has been done using both RM masses 



and SE masses (e.g., !Laodll998t iGreene & Holl2006t iBentz et al.ll2009cl: IXiao et al.ll2011 ). As 



suming some virial coefficient (/), these studies were able to add active objects in these scaling 
relations and extend the dynamic range in BH mass. 



In the past a few years, there have been a huge amount of effort to quantify the evolution of 
these scaling relations up to z ~ 6, by measuring host properties in broad-line quasars. Some 
studies directly measure the galaxy prop erties by de composing the quasar an d galaxy light in 



either imaging o r spectroscopic da t a (e.g..lTreu et al. 1|2004; Peng et aljbooeallbllWoo et aljbooel: 
Treuetal. 2007; Woo et al. 2008; Shen et al. 2008b; Jahnke et al.' 2009; 'McLeod & BechtoL 



2009:.Decarli et al...20 10; Merloni et al...2010-.Bennert et al.,.2010;.Cisternas et al..201 l:,Targett et al 



2OI2I) : other studies use indirect methods to infer galax y properties, such as using the narrow 
emission line width to infer bulge v elocity dispersion ('e.g.. lShields et al .I2OO3I l2006tlSalviander et al. 
2007t ISalviander & Shieldsll2012h . Molecular gas (using CO tracers) has also been detected in 
the hosts of z ~ 6 qua sars, allowing rough estimates on the host dynamical mass of these highest 
redshift quasars (e.g.. IWalter et al.ll2004t IWang et al.luOlOi and references therein). In all cases 
the BH masses were estimated using the SE methods based on different broad emission lines. 
With a few exceptions, most of these studies claim an excess in BH mass relative to bulge proper- 
ties either in the Mbh - cr* relation or in the Mbh - A/buige/i'buige relation, and advocate a scenario 
where BH growth precedes spheroid assembly. 
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It is worth noting that measuring host galaxy properties of Type 1 AGNs could be challen- 
ging, and systematic biases may arise when measuring the stellar velocity dispersion fr om s pectra 
(e.g.jBennert et al.ll201 lb . or host luminosities from image decomposing (e.g., iKim et al] 20081 : 
Simmons & UrrviEOOSl) . Conversions from measurables (such as host luminosity) to derived 
quantities (such as stellar mass) are also likely subject to systematics, especially for low-quality 
data. Thus careful treatments are required to derive unbiased host measurements. 

On the other hand, it is also worrisome that the errors in SE BH mass estimates may af- 
fect the observed offset in the BH scaling relations at high-redshift. As discussed extensively in 
^ there are both physical and practical concerns that the applications of locally-calibrated SE 
estimators to high-redshift quasars may cause systematic biases. Even if the extrapolations are 
valid, the luminosity-dependent bias discussed in 93. 3 .21 may still l ead to an average overestim- 
ation of quasar BH masses in flux-limited surveys. IShen & Kellvl (1201 Ol) studied the impact of 
the luminosity-dependent bias on flux-limited quasar samples, and found an "observed" BH mass 
excess of ~ 0.2 - 0.3 dex for Lboi S 10'*^ ergs"' with a reasonable value of crj - 0.4 dex (see 
93.3.2l for details). This sample bias using SE mass estimates becomes larger (smaller) at higher 
(lower) threshold quasar luminosities. 



Another statistical bias was pointed out by iLauer et al.l (120071) . which is at work even if 
there is no error in BH mass estimates. The basic idea is that since there is an intrinsic scatter 
between BH mass and bulge properties (~ 0.3 dex for the local sample), and since the distribution 
functions in BH mass and galaxy properties are expected to be bottom-heavy, a statistical excess 
(bias) in the average BH mass relative to bulge properties arises when the sample is selected based 
on BH mass (or based on quasar luminosity, assuming the Eddington ratio is const ant). This is 



similar to the Mahnquist- type bias discussed in 93.3.11 One can work out (e.g., iLauer et al 



2007UShen & Kellvll20inh that the BH mass off^set introduced by this bias depends on the slope 
of the galaxy distribution function, as well as the scatter in the BH-host scaling relations. For 
simple power-law models of the galaxy distribution function on property S with a slope y,, and 
lognormal scatter cr^ at fixed galaxy property S , this bias takes a similar form as the Malmquist- 
type bias in 93.3.11 



AlogMBH,Lauei- = - ln(10)r,cr;;/C , 



(21) 



where C is the coeflicient of the mean BH-host property (S) scaling relation log Mbh - C log S + 
C . Thus if the intrinsic scatter in the BH-host scaling relations increases with redshift, then 
this statistical bias alone can co ntribute a sign ificant amount to the observed BH mass offset in 



the high-redshift samples (e.g.. iMerloni et al.ii2010.) . A larger intrinsic scatter in these scaling 



relations at high redshift is expected, if the tightness of the local BH-host scaling relations is 
mainly estab lished via the hierarchical merging of less-correlated BH-ho st systems at higher 
redshift (e.g., |Pengll2007l: iHirschmann et al ] |2010H iJahnke & MaccidEoT 1|). The real situation is 
of course more complicated, and one must consider a realistic Eddington ratio distribution at fixed 
BH mass and the effect of the flux limit. There could also be other factors that may complicate 
the usa ge of AGNs to probe the ev olution of these BH-host scaling relations, as discussed in 
detail in lSchulze & Wisotzkil ( 201 1 ). But overall a BH mass excess due to the Lauer et al. bias is 
expected when select on quasar luminosity. An interesting corollary is that a deficit in BH mass 
is expected if the sample is selected based on galaxy properties. This may explain the findings 
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that high-redshift submillimeter galaxies (SMGs) tend to have on average smaller BHs relative to 
expectations from local BH-host scaling relations (e.g., Alexander et al. 2008). 

The Lauer et al. bias caused by the intrinsic scatter in BH-host scaling relations works in- 
dependently with the luminosity-dependent bias caused by errors in SE masses, so together they 
can c ontribute a sub stantial (or even full) amount of the observed BH mass excess at high redshift 
(e.g., Shen & Kellvl[201Q) . Both biases are generally worse for sarrmles with a higher luminosity 
threshold given the curvature in the underlying distribution function^ thus higher-z samples with 
higher intrinsic AGN luminosities will have larger BH mass biases, leading to an apparent evolu- 
tion. There are s everal samples that are probing similar luminosities as the local RM AGN sample 
(e.g.. lWoo et"aDl2 006). Since the SE mass estimators were calibrated on the local Mbh - cr* rela- 
tion using the RM AGN sample, one argument often made is that both biases should be calibrated 
away for the high-z sample with similar AGN luminosities. This argument is flawed, however, 
because the local RM AGN sample is heterogeneous and is not sampling uniformly from the 
underlying BH/galaxy distribution functions, while the high-z sample usually is sampling uni- 
formly from the underlying distributions - this is exactly why both biases will arise for the high-z 
samples. The only exception that might work is to compare two quasar samples at two different 
redshifts with the same luminosity threshold, where the predicted BH mass biases should be of 
the same amount, and see if there is evolution in the average host properties. But even in this 
case it requires that the underlying distributions (slope and scatter) and measurement systematics 
are the same for both the low-z and high-z samples. Proper simulations that take into account 
the measurement systematics (in both BH mass and host properties) and underlying distributions 
should be performed to verify the interpretations upon the observations. 

To summarize, there might be true evolutions in the BH-host scaling relation^^ but the 
current observations are inconclusive, due to the unknown systematics in the BH mass and host 
galaxy measurements. Better understandings of these systematics, the selection effects, as well as 
theoretical priors are all ne eded to probe the e volution of these scaling relations, and such efforts 
have been underway (e.g., CrotonI 2006: Rob ertson et al. 2006 : Lauer et al. 2007 : Hopkins et al 



2007; Di Matteo et al. 2008; Booth & Schave 2009"; Shankar et al. 2009; Shen & Keflv '201^ ; 
Hirschmann et al..2010; Schulze & Wisotzki201 1; Jahnke & Maccia201 l;,Portinari et al...2012: 



Salviander & Shieldsll20lilZhang et al.ll2012hT 



5. Summary and Future Perspectives 

To conclude, there have been considerable progress over the past several decades on the de- 
velopment of BH weighing methods for quasars. We now have a working technique based on 



This is not always tlie case. If tlie scatter (cr^, crj) increases at the low-mass end, then both biases could be worse at 
the low-mass/luminosity end. 

'''if the mass and velocity dispersion of galaxies bear any resemblance to the virial mass (Mh vir) and virial velocity 
(^^h.vir) of their host dai'k matter halos, then one of the two relations, A/bh ~ o"» and Mbh ~ !^bu\gc, must evolve since the 
^h.vir ~ ^h.vii- relation is redshift-dependent. 
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reverberation mapping of broad line AGNs that can measure active virial BH masses with an 
accuracy of a factor of a few (~ 0.5 dex). Rooted in the RM technique, efficient SE virial mass 
estimators have been developed to measure BH mass for large statistical samples of broad line 
quasars based on single-epoch spectroscopy. These methods greatly facilitate quasar studies in 
the era of modem, large-scale spectroscopic surveys. 

However, there are outstanding issues regarding the reliability of these virial (RM or SE) 
mass estimators, and consequences of their significant uncertainties, which is the focus of this 
review. Specifically I have the following highlighting remarks. 



• There are genuine concerns that the current RM sample does not represent the whole 
quasar/ AGN population, and the limited sample size and luminosity/redshift ranges of the 
RM sample, as well as the poorly understood BLR structure and dynamics, may impact 
the applicability of these locally-calibrated SE relations to high-z and/or high-luminosity 
quasars 03. lb . Due to limitations of the current RM sample and the uncertainty in the 
average virial coefficient (/), systematic biases on the order of a factor of a few are likely 
present due to these physical caveats. 



Even when the extrapolation of these SE estimators to high-z quasars is justified, rigor- 
ous statistical biases will arise from the uncertainties (scatter) in these SE masses ( ^3. 3b . 
In particular I demonstrated the conceptual difference between errors in SE masses and 
the distribution of SE masses within restricted luminosity ranges. Since luminosity is 
used in the estimation of SE masses, the variance in SE mass is redu ced when luminos- 
ity is constrained. I also derived the luminosity-dependent bias (e.g.. IShen et al.ll2008a : 



Shen & Kellvl201ol 2012 ) that at fixed true mass, the SE masses are over(under)-estimated 



in the mean when luminosity is higher (lower) than the mean luminosity at this fixed true 
mass, due to the stochastic variations between luminosity and line width that contribute to 
the uncertainty of SE masses. Simple simulations were performed to demonstrate these 
effects, and suggest that sample biases on the order of a factor of a few are present in 
flux-limited bright quasar sample. Thus these error-induced biases are as significant as the 
unknown systematic biases in SE masses, and cannot be ignored. More importantly, even 
when we eliminate all systematic biases (zero-point uncertainty) of these SE estimators 
in the future, these error-induced statistical biases will largely remain given the imperfect 
nature of these estimators. The formalism in ^3.3.2l provides guidance on how to quantify 
these error-induced sample biases with simulations. 



• Properly accounting for the selection effect due to the sample flux limit, and statistical bi- 
ases arising from errors in SE masses, are crucial to interpreting the observed distributions 
of quasars in statistical samples (@. I discussed how the "observed" distributions of BH 
mass and Eddington ratio for threshold data and with SE masses differ from the "true" 
distributions ( ^4. 2b . and cautioned on some recent claims based directly on the observed 
distributions. I further commented on the impact of the error-induced biases in SE masses 
on studies of the evolution of the BH-host scaling relations ( ^4. 3b . and concluded that the 
current observations are inconclusive for the claimed evolution. 
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Looking forward, there is an urgent need to expand the current RM sample with lag meas- 
urements, and, to acquire exquisite (velocity-resolved) RM data to utilize the full power of this 
technique. Only with a substantially larger RM sample that properly probes the AGN parameter 
space, and with a much better understanding of the BLR geometry and dynamics (for different 
lines) based on these RM data and ancillary data, can we improve these virial BH mass estimators 
further Since resource-wise, RM is a consuming exercise, it would be interesting to explore the 
possibilities of more efficient RM with wide-field multi-object imaging and spectroscopy, as well 
as dedicated facilities for single-object-mode monitoring. 



There are a few recent innovative proposals regarding RM that are worth mentioning. IZu et al, 



(1201 11) proposed an alternative method to measure RM lags, by fitting the observed continuum 
and emission line light curves with recently-developed statis tical models to describe quasar vari- 
ability (e.g., the "damped random walk" model developed bv lKellv et al ]|2009bl: iKozlowski et alj 
2010l). Compared with t he traditional cross-correlation method in measuring RM lags (e.g.. 



Gaskell & Petersonlll987h . this new method can improve lag measurements by simultaneously 
fitting multiple lines and quantifying error correlations. In addition, there have been efforts 
to bui ld dynamical BLR models to directly m odel the RM light curves in the time domain 



(e.g.. iBrewer et al 



201 It iPancoast et al.ll2012l) : such forward modeling (as discussed in, e.g. 



Netzer & Petersonlll997h can in principle provide direct constraints on the geometry of the BLR 
and the mass of the BH, and so it is worthwhile to explore its potential fur ther. Finally, al- 



Haas et al. 2011 



ternative strategies in RM experiments with no or few spectroscopic data (e.g 

Chelouche & DanieJ2012H Fine et al ]l2012 h. while not as good and reliable as traditional spectro- 
scopic RM, may speed up the process of probing the diversity of AGN parameters in the context 
ofRM. 
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